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SECTION  1. 
INTRODUCTION 


1.1.  Background 

The  design,  analysis,  and  fabrication  of  metal  naval  structures  has  matured  to  the  point 
that  fairly  explicit  design  rules,  factors  of  safety,  and  failure  prediction  methodology  exist 
in  the  form  of  handbooks  and  design  data  sheets.  Application  of  fiber  reinforced 
composite  materials  to  naval  structures  is  considered  a  “new”  technology,  although  the 
aircraft  and  pleasure  (small)  boat  industry  have  used  composite  structural  components  for 
more  than  thirty  years.  Composite  technologies  from  these  industries  have  been  slow  to 
transition  to  naval  applications.  Composite  small  boat  structures  are  typically  loaded 
differently  and  far  less  severely  than  naval  warships.  Aircraft  applications  typically  use 
composite  laminates  much  thinner  than  those  required  by  naval  vessels.  Furthermore, 
aircraft  fabrication  techniques  for  composites  are  economically  unfeasible  for  large  naval 
components  such  as  hulls.  Decision-makers  within  the  Navy  are  uncomfortable  with  a 
general  lack  of  and  inconsistency  in  data  associated  with  the  design,  analysis,  and 
fabrication  of  marine  composite  structures. 

The  Navy  has  attempted  to  resolve  technical  uncertainties  considered  unique  to  marine 
applications  of  composites.  New  and  much  improved  fabrication  methods  have  become 
available  [1]  reducing  the  variability  of  mechanical  properties  in  marine  composites. 
Commercial  finite  element  programs  have  introduced  new  elements  and  analysis 
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techniques  aimed  specifically  at  composite  structures,  facilitating  more  detailed  and 
accurate  structural  response  predictions.  As  these  technologies  advanced,  naval 
composite  structural  designs  moved  from  concept  to  prototype  stages.  Once  these 
composite  applications  began  to  appear  on  naval  vessels,  the  capability  to  analytically 
predict  composite  structural  failure  became  a  major  issue. 

The  ever-growing  size  of  many  Navy  composite  projects  such  as  masts,  deckhouses, 
hangers  and  stacks,  dictates  that  associated  components  become  relatively  thick.  Use  of 
fiber-reinforced  plastic  lamina  as  a  basic  structural  material  results  in  creation  of 
laminates  having  ten  or  more  lamina.  To  make  matters  more  complex,  each  lamina  may 
have  a  different  orientation  to  tailor  a  laminate’s  strength,  stiffness  and  functionality  to 
load  and  service  requirements.  To  analyze  new  designs,  engineers  compensate  for  this 
increased  complexity  by  smearing  (combining)  the  constituents’  (fiber  and  matnx) 
mechanical  properties  into  an  assumed  homogeneous  lamina.  They  then  smear  the 
lamina’s  mechanical  properties  into  an  assumed  homogeneous  laminate  to  attain 
idealized  uniform  strength  and  stiffness  throughout  the  structure.  While  this  approach  is 
fairly  successful  in  developing  global  stiffness  properties,  it  masks  large  internal  stresses 
associated  with  fiber  and  matrix  interactions.  Thus  the  designer  cannot  accurately  predict 
pre-,  ongoing,  and  post-damage  conditions  of  the  laminate  at  the  constituent  level  where 
damage  occurs.  At  the  other  extreme,  researchers  have  attempted  to  integrate  complex 
micromechanics  models  into  structural  analysis.  The  iterative  nature  of  this  solution 
process  results  in  excessive  computational  time  and  data. 
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Predicting  failure  of  structural  composites  has  been  a  topic  of  research  since  the  mid- 
1 960’s.  A  review  of  some  of  that  research  pertinent  to  this  dissertation  is  presented  later. 
My  experience  has  been  that  very  few  failure  prediction  methodologies  developed  in 
academia  have  been  used  in  Navy  projects  developed  by  government  and  private 
industry.  There  are  many  possible  explanations  why  the  methods  failed  to  make  the 
transition:  they  were  perceived  as  too  complicated  in  view  of  the  amount  of  time  and 
funding  allocated  to  the  problem;  they  were  formulated  around  plane  stress  assumptions 
and  could  not  be  applied  to  general  three  dimensional  stress  states;  or  thick  laminates 
(typical  in  navy  structures)  simply  presented  too  big  a  numerical  problem. 

Typically,  failure  analyses  assume  a  homogeneous  laminate  and  use  the  maximum 
attainable  stress  or  strain  values  of  the  laminate  as  determined  from  uniaxial  laboratory 
tests  as  critical  failure  parameters.  Among  the  shortcomings  of  this  method  are  that 
mitigating  or  exacerbating  effects  of  interactions  among  the  six  stress  components  are 
ignored,  all  possible  uniaxial  stress  states  are  not  usually  tested  (typically  through 
thickness  strengths),  and  strength  values  are  unique  to  the  architecture  of  a  particular 
laminate.  If  any  of  the  lamina  orientations,  groupings,  or  thicknesses  are  changed, 
maximum  laminate  stress  and  strain  values  change  also  and  new  tests  must  be  conducted. 

When  compared  to  the  analysis  of  metallic  structures,  lack  of  a  general  failure  prediction 
methodology  for  composites  leads  to  mixed  analytical  results  and  low  confidence  failure 
predictions.  Structural  designers  and  program  managers  use  large  factors-of-safety  to 
offset  a  lack  of  prediction  accuracy.  A  result  of  this  strategy  is  that  composite  structures 
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are  heavier  and  cost  more  than  necessary.  In  an  environment  in  which  composite 
structures  are  viewed  as  alternates  or  competitors  to  metal  structures,  over-design  causes 
them  to  compare  less  favorably  and  hinders  their  introduction  into  everyday  use. 

1.2.  Research  Objective 

The  intent  of  this  research  effort  is  to  develop  a  design  tool  for  failure  analysis  of 
composite  structural  laminates  that: 

a)  Conforms  to  analysis  methods  already  in  use, 

b)  Is  computationally  fast, 

c)  Produces  accurate  results,  and 

d)  Requires  a  minimum  of  additional  work,  on  the  part  of  the  user,  to  prepare  analysis 
inputs. 

These  requirements  tend  to  be  contradictory.  Computational  speed  is  often  accomplished 
by  simplifying  assumptions  that  reduce  accuracy.  Restricting  analysis  methods  to  current 
practice  and  minimizing  data  input  also  tend  to  oppose  a  goal  of  high  accuracy.  Hence, 
creating  a  design  tool  that  balances  these  requirements  constantly  requires  compromising 
one  or  the  other  in  an  attempt  to  arrive  at  some  optimum  solution.  Such  an  approach  will 
always  be  open  to  questions  about  the  legitimacy  of  the  methodology. 

Recent  introduction  of  Multicontinuum  Theory  (MCT),  developed  at  the  University  of 
Wyoming,  offers  an  approach  to  potentially  meet  these  requirements.  MCT  is  a  theory 
and  associated  numerical  algorithm  for  extracting,  virtually  without  a  time  penalty,  the 
stress  and  strain  fields  for  a  composite’s  constituents  during  a  routine  structural  finite 
element  analysis  (FEA).  Fully  implemented  in  a  FEA  program,  MCT  provides  the 
designer  with  additional  insight  into  a  composite  structure’s  response  to  its  load. 
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c 

Constituent  information  is  valuable  to  analysts  interested  in  predicting  composites’ 
failure.  For  instance,  advanced  failure  theories  for  composites  may  be  implemented  in  a 
more  straightforward  and  general  manner  at  a  constituent  level.  Damage  in  a  composite 
typically  begins  at  the  constituent  level  and  may  in  fact  be  limited  to  only  one  constituent 
in  some  situations.  An  accurate  prediction  of  constituent  damage  at  a  point  in  a  laminate 
provides  a  genesis  for  progressively  analyzing  failure  of  composite  structure  from  start  to 
finish. 

1.3.  Research  Overview 

The  research  contained  herein  concerns  the  development  of  a  constituent-based  failure 
criterion  and  its  incorporation  into  a  nonlinear  finite  element  code  to  produce  a  design 
tool  with  progressive  failure  prediction  capability.  Progressive  failure  analysis  in  an 
MCT-based  finite  element  code  has  opened  the  possibility  of  accurate  analysis  of 
composite  structures  in  pre-,  ongoing  and  post-  damage  conditions.  Because  information 
is  developed  on  a  level  at  which  failure  initiates,  load  redistribution  to  other  parts  of  the 
structure  as  well  as  the  remaining  constituents  can  be  efficiently  included.  This  allows  a 
designer  to  track  failure  as  it  occurs,  region  by  region,  and  reduce  the  stiffness  and 
strength  of  damaged  areas  without  necessarily  declaring  total  structural  failure.  Changes 
in  the  structure,  to  remedy  weaknesses,  can  be  concentrated  to  specific  areas  leading  to  a 
more  optimized  design. 
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Section  2  is  a  modest  review  of  the  published  literature  concerned  with  micromechanics- 
based  structural  failure  analysis.  Section  3  develops  the  equations  of  Multicontinuum 
Theory.  Section  4  develops  the  constituent-based  failure  criteria  used  within  this 
research.  Section  5  discusses  the  importance  in  developing  tightly  bound  constituent  and 
composite  properties  and  the  often  neglected  influence  of  nonlinear  shear  response. 
Section  6  describes  the  MCT_77  finite  element  code  and  implementation  of  the 
mechanics  issues  discussed  in  Sections  three  through  five.  In  Section  7,  comparisons  of 
analytical  results  against  experimental  data  are  provided  to  assess  the  validity  of  the 
methodology,  which  is  summarized  in  Section  8  discussions.  Stress-strain  graphs  and 
tables  of  important  parameters  for  the  materials  used  in  this  research  are  found  in  the 
Appendices. 


6 


NSWCCD-65-TR-1999/15 


SECTION  2. 

LITERATURE  REVIEW 


Most  failure  criteria  developed  for  composite  structural  laminates  to  date  can  be 
classified  as  macromechanical  because  the  criteria  attempt  to  predict  failure  using 
composite  stress-strain  data.  A  key  element  of  macromechanics  is  the  combining  of  the 
constituent’s  properties  into  a  homogeneous  set  of  composite  lamina  properties  and 
possibly  combining  lamina  properties  in  homogeneous  laminate  properties.  Reviews  of 
macromechanical  failure  criteria  can  be  found  in  Tsai  [2],  Rowlands  [3],  and  Nahas  [4]. 

In  contrast,  micromechanical  approaches  retain  the  individual  identities  of  each  lamina 
and  its  constituents.  This  allows  interaction  among  them.  Composite  properties  are  used 
in  micromechanics  analyses  but  failure  of  each  constituent  and  its  contribution  to  lamina 
and  laminate  failure  is  emphasized.  All  micromechanical  models  are  predicated  on  a 
complete  set  of  material  constants  for  each  constituent.  Composite  properties  for  a 
unidirectional  lamina  are  typically  synthesized  from  constituent  properties  via  a  finite 
element  or  closed  form  analytical  model. 

MCT  is  micromechanics  based  and  routinely  generates  constituent  level  stress  and  strain 
fields  during  an  analysis.  Hence,  macromechanical  based  failure  criteria  and  analyses  are 
not  of  primary  interest  in  the  present  research.  Instead,  the  focus  is  on  constituent  based 
failure  approaches.  Such  approaches  may  provide  higher  resolution  in  predicting 
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composite  material  failures  since  failure  often  begins  in  one  constituent  (e.g.  matrix 
cracking)  while  the  other  constituent  remains  intact. 


Analysis  of  structural  laminate  failure  based  on  micromechanics  has  been  largely  ignored 
in  the  past  because  constituent  information  was  generally  unavailable  in  standard  finite 
element  structural  analysis.  There  were,  however,  notable  exceptions.  Aboudi  [5] 
extended  a  simplified  micromechanical  model  he  developed  earlier  to  predict  failure 
surfaces  for  unidirectional  fibrous  composites.  The  micromechanics  model,  known  as  the 
“Method  of  Cells”,  was  based  on  a  representative  unit  “cell”  of  a  composite  assumed  to 
be  composed  of  a  periodic  array  of  square  fibers  embedded  in  a  matrix  (Fig.  2.1  and  Fig. 
2.2).  The  model  allowed  for  linear  variation  of  displacements  within  each  of  four 
subcells  making  up  a  representative  unit  cell.  Continuity  of  average  displacements  and 
tractions  was  enforced  between  subcells.  Relationships  between  composite  stresses  and 
those  of  the  constituent  subcells  were  obtained  in  closed- form.  This  closed-form 
relationship  allowed  for  efficient  constituent  stress  calculations  from  structural  level 
results,  thus  facilitating  implementation  of  micromechanics-based  failure  criteria. 


Fig.  2.1  Aboudi's  idealized  unidirectional  fiber  composite. 
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Fig.  2.2  Aboudi's  representative  unit  cell. 


Three  types  of  failures  were  assumed:  fiber,  matrix,  and  fiber/matrix  interface.  It  was 
noted  that  constituent  in  situ  ultimate  strengths  were  generally  different  from  bulk 
measured  values  and  that  values  for  constituent  strengths  are  different  in  compression 
versus  tension.  Aboudi  theorized  that  if  ultimate  strengths  of  a  unidirectional  composite 
in  principal  material  directions  could  be  determined,  ultimate  in  situ  constituent  strengths 
could  be  determined  by  backing  them  out  via  the  micromechanical  model.  Failure 
envelopes  could  then  be  generated  for  a  composite  under  any  complex  stress  state. 
Aboudi  used  this  approach  to  implement  an  uncoupled  (fiber  from  matrix)  maximum 
stress  criterion  as  follows.  Retaining  a  plane  stress  assumption  used  at  the  lamina  level, 
composite  failure  occurred  when  one  of  the  following  conditions  was  satisfied: 


22 


>+s * 


12 


>s: 


where: 

+Si/=  Fiber  tensile  strength, 
Ssm  =  Matrix  shear  strength, 

022 m  =  Transverse  matrix  stress, 


+S  m  =  Matrix  tensile  strength, 
oi/~  Axial  fiber  stress, 

0/2  m  -  Matrix  in-plane  shear  stress 
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The  damage  algorithm  implied  by  the  above  criterion  is  one  of  catastrophic  failure  in 
both  the  constituent  and  composite  when  either  of  the  constituents  is  damaged.  This  type 
of  damage  algorithm  is  perhaps  overly  simplistic  in  light  of  experimental  evidence 
showing  that  matrix  cracking  is  not  necessarily  catastrophic  for  the  composite. 

To  verify  the  micromechanics  model  failure  prediction  capabilities,  Aboudi  made 
extensive  comparisons  of  predicted  off-axis  tensile  strengths  of  various  unidirectional 
lamina  with  experimental  data  from  six  composites:  boron/epoxy,  graphite/polyimide, 
AS4/3501  graphite/epoxy,  aramid/epoxy,  glass/epoxy,  and  borOn/aluminum.  There  was 
good  agreement  between  predicted  and  measured  strengths.  Fig.  2.3  shows  a  typical 
result.  Failures  in  these  composites  were  catastrophic  when  either  constituent  was 
damaged  so  Aboudi’s  failure  criterion  worked  well  in  these  cases. 

Pecknold  [6]  and  Rahman  [7]  extended  a  micromechanical  model  developed  earlier  by 
Pecknold  to  predict  composite  structural  failure  within  the  framework  of  a  finite  element 
program.  Their  micromechanical  model,  like  that  of  Aboudi,  assumes  a  composite 
microstructure  composed  of  square  fibers  arranged  in  a  doubly  periodic  array  (Fig.  2.1). 
The  model,  referred  to  as  a  “Combination  Type”,  was  based  on  a  mechanics  of  materials 
approach  and  recognized  the  restraining  effect  of  the  fibers  on  the  matrix.  The  model 
used  a  representative  unit  cell  similar  to  that  used  by  Aboudi  (Fig.  2.2),  except  two  matrix 
subcells  were  combined  for  a  total  of  three  instead  of  four  subcells.  Effective  elastic 
constants  for  the  composite  unit  cell  were  developed  by  use  of  spring  analogies 
representing  a  subcells’  stiffness.  Constant  stress  or  constant  strain  was  assumed  within 
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each  constituent,  depending  on  load  direction.  This  approach  led  to  a  simpler 
formulation  of  composite  properties  than  the  Method  of  Cells,  but  with  comparable 
accuracy.  The  micromechanics  model  was  linked  to  and  drove  the  composite  material 
response  of  the  finite  element  analysis  at  each  material  (Gauss)  point. 


Fig.  2.3  Aboudi ’s  method-of  cells  comparison  between  predicted  and  measured  off-axis 
strength  of  a  glass/epoxy  composite  [5], 

Rahman  recognized  six  failure  criteria:  five  uncoupled,  in-plane  only,  maximum  stress 
micromechanical  criteria  and  one  macromechanical  criterion.  Failure  occurred  when  one 
of  the  conditions  listed  in  Table  2.1  was  satisfied.  Implied  in  the  criteria  was 
directionality  (orientation)  of  the  failure,  e.g.,  matrix  cracks  along  the  fiber  axis  when 
transverse  tensile  failure  occurred. 
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Table  2.1  Failure  criteria  used  in  the  Combination  Type  micromechanical  model. 


Failure  Type 

Condition 

Fiber  fracture 

Matrix  shearing 

<2  *  W 

Transverse  matrix  cracking  in  tension 

+  rtm 

ct22>  b22 

Transverse  matrix  crushing 

-  m  -r»m 

<t22  —  o22 

Axial  matrix  cracking  in  tension 

+<>+S” 

Lamina  kink  banding  (fiber  micro-buckling) 

-cr„  >  (FVF)Gn 

In  Table  2.1: 

.  FVF  =  Composite  fiber  volume  fraction, 

G12  -  Composite  in-plane  shear  modulus, 

'S22  m  =  Matrix  in  situ  compression  strength, 

'(711  =  Compressive  stress  in  the  lamina, 

+<7nm  =  Axial  tensile  stress  in  the  matrix. 

Once  constituent  damage  was  detected,  constituent  material  stiffnesses  were  gradually 
degraded  and  relevant  stress  components  were  either  fully  relieved  (fiber  fracture  or 
matrix  cracking)  or  held  constant  (matrix  shearing,  crushing  or  lamina  kink  banding)  at 
the  Gauss  point  being  evaluated.  The  material  softening  and  stress  release  allowed  the 
load  to  be  redistributed  within  the  composite  to  surrounding  undamaged  areas.  A 
stiffness  degradation  algorithm  uncoupled  failure  direction  response  from  other  stress  and 
strain  components  within  the  compliance  matrix.  This  was  accomplished  by  setting  off- 
diagonal  terms  to  zero  and  increasing  the  diagonal  term  compliance  by  the  inverse  of  a 
specified  stiffness  reduction  factor.  As  an  example.  Fig.  2.4  shows  a  fiber  compliance 
matrix  after  axial  failure. 

The  stiffness  reduction  factor,  A>  controlled  how  rapidly  the  load  redistribution  occurred 
and  its  value,  empirically  determined,  was  chosen  to  improve  solution  convergence. 
Rahman  [7]  noted  that  solutions  appeared  insensitive  to  a  wide  range  of  stiffness 
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reduction  factors  and  included  examples  such  as  0.25  for  shear  modulus  and  0.01  to 
0.001  for  Young’s  modulus. 

Z)"'*,,  0  0  0  0  0 

S22  S23  S2A  S25  S26 

S33  ^34  ^35  S36 

SA4  ^45  *^46 

Sym  s55  sS6 

*66 

Fig.  2.4  Rahman 's  fiber  compliance  matrix  after  tensile  failure  [7J. 

Rahman  used  the  Combination  Type  micromechanics  model  and  associated  failure 
criteria  to  predict  and  compare  off-axis  lamina  strengths  of  the  same  experimental  data 
sets  used  by  Aboudi.  Results  were  in  excellent  agreement  with  experimental  data  and 
virtually  identical  to  those  from  Method  of  Cells.  Rahman  went  on  to  analyze  more 
complex  laminates  and  incorporated  nonlinear  shear  modulus  behavior  in  the 
micromechanical  model.  Excellent  agreement  between  analytical  results  and 
experimental  data  of  composite  failure  strengths  was  achieved  in  four  tensile-loaded, 
boron/epoxy  laminate  configurations  ([0/90]s,  [±30]s,  [±60]s,  [±45/903]s).  Good 
agreement  in  composite  failure  strength  and  correct  identification  of  failure  type  was 
achieved  in  the  analysis  of  four  compressive-loaded,  boron/epoxy  laminates  ([0/90]s, 
[±30]s,  [±60]s,  [±20]s). 

Rahman  used  the  failure  criteria  in  a  progressive  failure  algorithm  within  the  finite 
element  framework.  The  algorithm  was  used  to  simulate  the  compressive  behavior  of  a 
notched  (open  hole),  [±45]6s,  T300/BP976  laminate  to  failure.  The  algorithm  indicated 
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that  initial  damage  was  fiber  micro-buckling  around  the  hole’s  edge,  which  precipitated 
local  matrix  shear  failure.  A  reduction  in  shear  stiffness  followed  which  precipitated 
additional  fiber  micro-buckling  cascading  into  laminate  failure.  Analytical  results  were 
in  good  agreement  with  experimental  data  predicting  a  failure  load  11%  higher  than 
actual.  Rahman  noted  in  his  conclusions  that  a  possible  improvement  in  his  results, 
especially  in  predicting  compressive  failure,  may  have  been  achieved  by  using  interactive 
failure  criteria  that  account  for  mitigating  and  exacerbating  reactions  between  the  stresses 
rather  than  simple  maximum-stress  criteria. 

Kwon  [8]  also  developed  a  micromechanics  model  to  describe  failure  and  damage 
progress  in  laminated  fibrous  composite  structures.  Kwon’s  model  assumed  the  same 
periodic  composite  microstructure  and  unit  cell  geometry  as  Aboudi  (Fig.  2.1  and  Fig. 
2.2).  Like  Aboudi,  constant  stress  and  strain  is  assumed  within  each  subcell,  which 
simplifies  determination  of  constituent  (subcell)  stresses  from  those  of  the  composite 
(closed-form  equations).  The  model  was  incorporated  into  a  finite  element  program  but, 
unlike  Pecknold  and  Rahman,  it  does  not  drive  the  composite  material  response  of  the 
structure.  Instead,  the  model  uses  composite  strains  computed  in  the  course  of  normal 
structural  analysis  (assuming  pseudo-homogenous  material  properties)  to  derive 
constituent  strains  and  stresses. 

Kwon’s  micromechanical  failure  criteria,  like  that  of  Aboudi,  Pecknold,  and  Rahman, 
was  based  on  fiber  and  matrix  stresses.  Assuming  plane  stress,  failure  occurred  when  one 
of  four  conditions  was  satisfied  in  a  local  material  coordinate  system.  For  fiber  breakage 
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f  =  1  if  <?{\  >  0  or 


if  a(x  <  0 . 


For  matrix  cracking  a  quadratic  interactive  failure  criterion  (the  form  noted  to  be  similar 
to  those  proposed  by  Hashin  and  Rotem,  [9]  and  Yamada  and  Sun  [10])  was  used: 


f  —.m 
gl2 

Sm 

s  J 


=  1 


if  <2  >  0 


or 


22 


Sm 

22  J 


'12 


1  if  cr”  <  0 


\as  J 


Once  constituent  failure  was  detected,  the  damage  algorithm  degraded  appropriate 
constituent  material  properties  by  setting  normal  and  shear  moduli  to  near  zero  values. 
Near  zero  values,  rather  than  zero,  were  necessary  to  avoid  numerical  difficulties 
associated  with  matrix  singularities. 


Kwon’s  analytical  predictions  of  off-axis  lamina  strength  were  in  excellent  agreement 
with  the  same  experimental  data  sets  used  by  Aboudi  and  Pecknold.  Kwon  extended  his 
damage  algorithm  to  analyze  progressive  damage  in  a  [(0/90)6]s,  graphite/epoxy,  cross- 
ply,  composite  laminate  containing  a  hole  under  a  tensile  load.  Analysis  was  in  excellent 
agreement  with  experimental  data  for  ultimate  failure  strength  (345  MPa  versus  360  MPa 
experimental)  as  well  as  type,  location  and  failure  path  direction  (66°-78.5°  analytical 
versus  75°  experimental  off  the  applied  load  axis). 


Micromechanics-based  analysis  of  fibrous  composites  is  not  limited  to  unidirectional 
lamina.  Blackketter  [11]  developed  a  finite  element  based,  three-dimensional 
micromechanics  model  for  describing  damage  propagation  in  a  plain  weave 


15 


NSWCCD-65-TR- 1 999/15 


graphite/epoxy  (AS4/3501)  composite  (Fig.  2.5).  Constituents  of  a  unit  cell  were 
identified  as  matrix  or  yarn.  A  separate  two-dimensional  finite  element  micromechanics 
model  (single  fiber  in  a  matrix)  was  used  to  generate  material  properties  of  yarn  from 
fiber  and  matrix  properties. 


Fig.  2.5  Plain  weave  fabric  idealization  showing  yarn  pattern. 


Like  Rahman,  Blackketter  used  stiffness  reduction  factors  to  control  load  redistribution  to 
undamaged  areas  and  to  account  for  both  the  type  and  direction  (orientation)  of  the 
failure.  Table  2.2  shows,  for  example,  a  transverse  tension  failure  would  result  in  a  near 
zero  value  for  the  transverse  modulus,  an  80%  reduction  to  the  in-plane  shear  moduli, 
while  the  remaining  moduli  keep  their  original  value. 


Table  2.2  Blackketter’ s  failure  modes  and  degradation  values  of  elastic  properties. 


Failure  Mode 

m 

mam 

Longitudinal  tension  (+Sn) 

WSEM 

Transverse  tension  (+S 22) 

1.0 

IiMb 

0.2 

Transverse  tension  (*833) 

1.0 

ES 

1.0 

0.2 

Transverse  shear  (S23) 

wsm 

MB 

Longitudinal  shear  (S13) 

umm 

1.0 

0.01 

Longitudinal  shear  (S12) 

1.0 

0.01 

1.0 

1.0 

1.0 

0.01 
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A  significant  conclusion  appearing  in  Blackketter’s  work  noted  that  nonlinear  shear 
stress-strain  behavior  was  caused  principally  by  damage  propagation  rather  than 
nonlinear  matrix  material  behavior. 

The  research  cited  above  represents  significant  achievements  in  constituent-based 
structural  analysis  but  further  improvements  can  be  made.  Driving  composite  response 
from  a  constituent  level,  as  in  the  research  of  Rahman  and  Pecknold,  must  be  done 
efficiently  or  computational  time  penalties  could  become  significant.  Extensions  of  a 
finite  element  failure  code  to  other  constitutive  models,  such  as  viscoelasticity,  could 
impose  severe  computational  penalties  that  would  render  the  method  a  less  desirable 
analysis  alternative  to  mainstream  users.  Rahman  and  Pecknold  note  that  using  a  stress- 
interactive  failure  criterion  rather  than  a  simple  maximum  stress  (local  material 
coordinate  system)  criterion  may  further  improve  failure  predictions.  But  using  a  stress- 
interactive  constituent  failure  criterion,  assuming  plane  stress  in  the  constituents  such  as 
that  used  by  Kwon,  does  not  address  the  general  load  case.  Although  a  composite 
laminate  or  lamina  may  be  in  a  state  of  plane  stress,  its  constituents  are  generally  in  a 
fully  three-dimensional  stress  state. 

Approaches  such  as  Method  of  Cells  are  restricted  by  geometrical  considerations.  For 
example,  the  continuous  fiber  composite  model  assumes  a  square  packing  array  with 
square  fibers.  Brockenbrough  et  al.  [12]  have  shown  that  packing  arrangements  have  a 
pronounced  effect  on  inelastic  behavior  in  composites.  One  of  virtues  of  MCT  is  a  model 
geometry  may  be  developed  for  any  microstructure,  including  random  distributions  of 
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reinforcement.  This  capability  may  be  critical  when  attentions  are  turned  to  analyses  of 
other  composite  material  architectures  such  as  woven  fabrics,  particulate,  and  short 


(chopped)  fiber. 


Interestingly,  in  the  literature  reviewed,  no  attempts  were  made  to  compare  the  relative 
accuracy  between  micromechanics-based  failure  predictions  and  macromechanics  failure 
predictions.  This  comparison  would  seem  important  in  assessing  the  net  benefits  gained 
from  the  additional  work  required  to  implement  micromechanics  based  failure  analysis 
methods. 
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SECTION  3. 

MULTICONTINUUM  THEORY 

3.1.  Continuum  Definition  of  Stress 

Multicontinuum  theory  (MCT)  development  is  presented  in  detail  in  papers  by  Gamich 
[13]  and  Gamich  and  Hansen  [14,15].  It  is  summarized  here  as  a  convenience  and  to 
emphasize  concepts  important  to  implementing  a  constituent-based  failure  analysis. 

MCT  begins  with  a  continuum  definition  of  stress  at  a  point.  The  concept  of  stress  in 
homogeneous  materials,  such  as  steel,  is  a  familiar  one  to  most  engineers.  Yet,  if  looked 
at  on  a  micro-scale,  Fig.  3.1,  one  sees  that  the  “homogeneous”  material  is  hardly 
homogeneous.  It  is  obvious  that  stresses  will  vary  significantly  from  point  to  point  across 
the  material  [16]. 


19 


NS  WCCD-65-TR- 1 999/ 1 5 


The  homogenized  value  used  to  characterize  the  stress  tensor  at  this  point  is  derived  by 
taking  a  volume  average  of  all  stresses  in  the  region: 


(3-1) 


where  D  is  the  region  representing  the  continuum  point.  The  concept  of  a 
multicontinuum  simply  extends  this  concept  to  reflect  coexisting  materials  within  a 
continuum  point.  Consider  a  composite  material  with  two  clearly  identifiable 
constituents  as  shown  in  Fig.  3.2. 


Fig.  3.2  Composite  lamina  as  a  multicontinuum. 
Using  equation  (3.1)  for  each  constituent 


and 


dV, 


(3.2) 


(3.3) 
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where 

D  =  Df{jDm.  (3.4) 

Combining  equations  (3.1),  (3.2),  and  (3.3)  leads  to 

Z  =  4fZf+tmZm>  (3-5) 

where  ^and  are  the  volume  fractions  of  fiber  and  matrix  respectively. 

Likewise,  for  strains  we  have: 


£=<t>f  £f+<t>n,£m-  (3-6) 

It  is  important  to  note  the  averaging  process  that  results  in  these  equations.  That  is,  we 
are  not  concerned  with  stress  and  strain  variations  through  individual  constituents  within 
D  but  only  with  their  average  values.  This  is  an  information  compromise  that  separates 
structural  analysis  from  micromechanical  analysis.  Accounting  for  stress  variations 
throughout  every  fiber  at  every  material  point  in  even  a  modest  structure  is  simply  not 
possible  or  desirable.  In  contrast,  providing  constituent  average  stress  and  strain  fields 
opens  a  new  and  manageable  information  window  on  a  composite  material’s  response  to 
a  load  as  compared  to  a  classic  single  continuum  approach. 


3.2.  Determining  Constituent  Stresses  and  Strains 

Finite  element  solutions  of  structural  problems  produce  stress  and  strain  fields  of  the 
composite  thereby  providing  two  of  the  four  unknowns  in  equations  (3.5)  and  (3.6).  To 
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provide  closure  of  the  equations,  constitutive  relations  are  required  for  the  composite  as 
well  as  the  constituents.  Here  we  change  from  direct  tensor  to  contracted  matrix  notation 
and  write 


^11  C12  C13  C14  C15  C16 

<T22 

^22  ^23  ^24  ^25  ^26 

^22 

^33 

c33  c34  c35  c36 

^33 

< 

>  rr: 

> 

<J,2 

C44  ^45  C46 

r12 

<7,3 

SYM  c55  c56 

Tn 

.Cr23  . 

C66  _ 

7 23. 

For  composite  and  constituents  there  follows: 

KM KKKMM)> 

KMcJ(kMK }), 


(3.7) 


(3.8) 

(3.9) 
(3.10) 


where  [C]  and  [Cp]  represents  material  stiffness  matrices  and  {£„} and  {£&>}  are  thermal 
strains.  Substituting  equations  (3.8),  (3.9),  and  (3.10)  into  (3.5)  and  using  (3.6)  gives 

Solving  for  {Sf}  gives 


^  (3.12) 

<h 
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where 

{eo}  =  0{a},  {ef0}=9{af},  = 

In  the  above,  a  represents  the  coefficient  of  thermal  expansion  and  G  is  the  relative 
temperature.  Equation  (3.12)  may  be  rewritten  as 

W  =  '  (313) 

where 

[A]  =  -^([C\-[C,\y  ([C]-[C.]), 

*! 

and 

{a}  =  ([C]-[c/])-|([C]{a}-«>/[c/]{a/)-«>JCj{aJ). 

Note  that  [A]  is  a  matrix  relating  mechanical  constituent  strains,  while  the  vector  {a} 
relates  constituent  thermal  strains.  Typically  [C/],  [ Cm ],  {a/},  and  {a„ ,},  are  known 
material  properties  of  the  constituents.  [C]  and  {a}  of  the  composite  are  developed  from 
micromechanical  modeling.  Hence,  [A]  and  {a}  are  known  a  priori  and  could  be 
calculated  independent  of  the  finite  element  program. 

Determining  constituent  strains  from  composite  strains  requires  substituting  equation 
(3.13),  into  (3.6)  and  solving  for  {sm}  as 

({*}-»{«})•  (3-i4> 
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where  [I]  is  the  identity  matrix.  Substituting  (3.14)  back  into  (3.6)  yields  an  expression 
for  {£/}  as 

({*)-<*>„  kb-  <315> 

4>f 

Equations  (3.14)  and  (3.15)  allow  constituent  strains  to  be  to  be  calculated  from 
composite  strains  at  any  point  in  the  finite  element  model.  Using  equations  (3.9)  and 
(3.10),  constituent  stresses  can  be  calculated. 


3.3.  MCT  and  Material  Failure 

Clearly,  the  composite  and  constituent  stresses  and  strains  developed  previously  are 
based  on  average  values  at  a  continuum  point.  In  the  finite  element  method,  numerical 
integration  samples  stress,  strain,  and  material  continuum  values  at  Gauss  quadrature  or 
“material”  points.  Implied  in  the  derivation  of  MCT  and  its  finite  element 
implementation  is  that  both  fiber  and  matrix  exist  at  every  point  in  the  composite.  The 
concept  of  average  values  at  a  continuum  point  applies  to  material  properties  also,  i.e., 
there  is  no  variation  across  that  point  in  either  constituent.  Therefore,  when  a  constituent 
material  fails  at  a  Gauss  (continuum)  point  it  fails  absolutely  without  directionality  and 
only  the  other  constituent  remains.  This  is  in  sharp  contrast  to  some  of  the  work 
discussed  previously  in  the  literature  review,  which  incorporate  failure  models  that 
attempt  to  include  directional  dependence. 
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SECTION  4. 

FAILURE  CRITERIA 

4.1.  Overview  of  Quadratic  Failure  Criteria 

A  general  overview  of  current  failure  criteria  is  given  in  Appendix  G.  Here  we  are  only 
concerned  with  failure  criteria  that  account  for  the  interaction  effect  of  the  various 
stresses  on  a  material’s  strength. 


The  Maximum  Distortional  Energy,  or  von  Mises,  criterion  is  the  most  widely  used 
criterion  for  predicting  yield  points  in  isotropic  metals.  This  criterion  has  also  been  used 
to  predict  brittle  failure  in  isotropic  materials.  The  isotropic  von  Mises  failure  criterion  is 
a  special  case  of  a  general  form  called  quadratic  interaction  criterion,  so  named  because  it 
includes  terms  to  account  for  interaction  between  the  stress  components.  Hill  [17] 
generalized  the  von  Mises  criterion  to  include  the  effects  of  orthotropic  behavior.  The 
Hill  criterion  is  given  by 


-o-2)2  +^2(^2 +  FA^  -o-i )2 

+  ^23  (^23  +  C32  f  +  ^13  +  CT31  )'  +  FU  (^.2  +  ^21  )’  =  U 

where 


1  11 

2  F\  —  — —  H — - — 

1  r»2  n2  o2 

^11  ^22  °33 


IF  = 
12 


1 


2  F2  =-L  +  -L-J L 

4  r»2  n2  o2 

022  O33  Ojj 

2F  =F 

ci  2  ’ 

°13 


2  F 
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In  the  above,  q,  are  stresses  referenced  to  the  principal  axes  of  orthotropy  (principal 
material  directions),  while  S0  represent  failure  strengths.  As  presented,  Hill’s  criterion 
assumes  no  difference  between  tensile  and  compressive  yield  strengths. 


A  generalized  quadratic  interaction  criterion  form,  suggested  by  Gol’denblat  and  Kopnov 


[18]  and  proposed  by  Tsai  and  Wu  [19],  is  given  as 


Ftc T,-  +Fij(Tiaj  =1, 


(4.2) 


where  F,  and  Fy  are  experimentally  determined  strength  tensors  and  contracted  tensor 
notation  is  used  (ij  =  1  to  6).  Hoffman  [20]  has  suggested  that  the  linear  terms,  Fi,  are 
necessary  to  account  for  differences  in  tensile  and  compressive  strengths,  whereas  Tsai 
and  Wu  state  that  they  are  necessary  to  account  for  internal  stresses.  Expanding  equation 
(4.2)  for  an  orthotropic  material,  using  symmetry  of  the  strength  tensors,  assuming  no 
interaction  between  normal  and  shear  stresses,  and  assuming  negative  and  positive  shears 
are  equal  gives 


■fj<7 1  +  F2<j2  +  F3cr3  +  Fncr2  +  F22<72  +  F22a]  +  F^&l 

+  F55<j  2  +F66ct26  +  2Fncricr2  +  2Fi2(Ti<J3  +  2F23(72<j3  —1 , 
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fb =  IF  I1  “ P(F' +  Fl  ' ' p!  (F" +  >]  ’ 

Fa  =-^[l-Q(F,  +F,)-Q,(FU  +FJ]  , 

F»  +F,)-R,(F2,  +  Fm)]. 

ZK 

The  terms  accounting  for  coupling  effects  between  normal  stresses,  Fn,  F13,  and  F23, 
require  experimental  determination  using  biaxial  tests.  For  example,  using  biaxial 
tension  failure  loads  P,  Q,  and  R,  <Ji=Gz=P,  gj=G3=Q  and  05=05=/?.  The  remaining 
coefficients  can  be  experimentally  determined  from  simple  uniaxial  or  shear  tests.  The 
Tsai-Wu  criterion  appears  to  be  the  most  widely  used  quadratic  stress-interaction 
criterion  for  predicting  macro  failure  in  composite  materials. 

Hashin  [21]  developed  a  three-dimensional,  stress  interactive,  failure  criterion  for 
unidirectional  lamina  that  recognized  two  distinct  and  uncoupled  failure  modes.  While 
the  failure  criterion  itself  is  based  on  composite  stresses,  it  constructs  a  piecewise 
continuous  failure  form  based  on  constituent  failure  modes. 

Hashin  states  that  most  failure  criteria  that  account  for  stress  interaction  are  quadratic  in 
stress  and  that  choice  of  a  quadratic  polynomial  is  based  solely  on  curve  fitting 
considerations.  Quadratic  criteria  are  third  in  a  hierarchy  of  potential  approximations  of 
failure  surfaces  (constant  or  maximum  stress,  linearly  varying  stress  and  quadratic).  A 
quadratic  approximation,  he  notes,  “is  the  simplest  presentation  which  can  fit  the  data 
reasonably  well  and  in  view  of  the  significant  scatter  of  failure  test  data  it  hardly  seems 
worthwhile  to  employ  cubic  or  higher  approximations.” 
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A  unidirectional,  fiber  reinforced  composite  is  assumed  to  be  transversely  isotropic.  A 
local  orthogonal  coordinate  system  is  defined  in  which  the  fiber  axis  serves  as  the 
principal,  x,,  material  direction,  and  x2t  x3  the  transverse  directions.  Since  material 
constants  assume  transverse  isotropy,  stresses  can  be  expressed  as  values  invariant  under 
rotations  about  the  principal  material  direction.  The  failure  state  of  either  constituent 
material  can  be  expressed  in  terms  of  these  transversely  isotropic  stress  invariants. 
Although  Hashin  derived  the  invariants,  Hansen  [22],  in  development  of  an  anisotropic 
flow  rule  for  plastic  behavior  in  composite  materials,  presented  a  different  form  (both 
forms  attributed  to  Spencer)  used  here.  The  five  transversely  isotropic  stress  invariants 
are: 


A  —  ’ 

1 2  ~  &21  +  °33  > 

h  =  Cr22  +  &33  +  2<T23  ’  ^4'4^ 

14  —  <7U  +  <713  , 

1 5  :=Cr22&l2  +  ^33^13  +  ^12^13^23  * 

Hashin’s  choice  of  a  quadratic  form  eliminates  I5  from  appearing  in  the  failure  criteria. 
Therefore  the  most  general  form  for  a  quadratic  criterion  is 

AJi  +BJ?  +  A2I2  +  B2I\  +  Cl2IJ2  +  4/3  +  A4Ia  =  1  •  (4-5) 

Feng  [23],  in  an  analogous  manner,  developed  a  quadratic  interaction  failure  criteria 
based  on  transversely  isotropic  strain  invariants. 
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At  this  point,  it  is  instructive  to  compare  the  criterion  of  Tsai  and  Wu  with  that  of  Hashin. 
Tsai  and  Wu,  in  reference  [19]  presented  a  form  of  equation  (4.3)  for  transversely 
isotropic  composites  as 


Fi°\\ 


+  F2^22 


+  Cr33 


)  F\ I 1  +  F22  (^22  "*"Cr33  +  2<J23  )+  F44  (<J]2  +C7]3) 
^F12  (*^11^22  ^1 1  *-^*33  )  ^F23  {^22^33  ^23  )  —  ^  ' 


(4.6) 


Rewriting  equation  (4.6)  in  terms  of  the  transversely  isotropic  stress  invariants  gives 

F,I,  +F1I2+Fui;  +F11Ii+FJ,+2FnIlI1+2F11(l22  -I,)=  1, 
or, 

FJ\  +FJl  +F2*2  +2F23/22  +  2FI2/,/2  +(F22  -2 F23)/3  +F4I4  =1.  (4.7) 

Comparing  equation  (4.7)  to  equation  (4.5),  shows  that  the  Tsai-Wu  criterion  for 
transversely  isotropic  materials  and  the  Hashin  failure  criterion  have  the  same  functional 
form.  Their  difference  is  in  defining  the  coefficients  of  the  stress  terms.  The  Tsai-Wu 
equation  is  used  to  define  a  smooth  and  continuous  failure  surface  in  both  the  tension  and 
compression  regions  of  space.  As  a  result,  the  stress  tensors  are  functions  of  both  tensile 
and  compressive  composite  strengths.  Hashin  identified  two  composite  failure  modes; 
fiber  versus  matrix  influenced,  and  develops  separate  equations  based  on  the  failure  mode 
to  determine  a  failure  state.  Hashin  further  recognized  that  a  composite  typically  has 
different  ultimate  strengths  in  tension  and  compression,  so  both  fiber  and  matrix  failure 
criteria  have  tensile  and  compressive  subforms.  Hence,  the  coefficients  of  the  stress 
terms  are  functions  of  only  tension  or  compression  strengths  resulting  in  a  continuous  but 
not  smooth  stress  space  failure  surface. 
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4.2.  MCT  Failure  Criterion 

In  a  major  departure  from  Hashin’s  failure  criterion,  a  failure  criterion  was  developed  for 
each  constituent  rather  than  the  composite.  The  criterion  is  in  the  form  of  equation  (4.5) 
and  uses  constituent  stress  information  produced  by  MCT.  As  a  consequence,  the 
transversely  isotropic  stress  invariants,  defined  in  equation  (4.4),  were  used  for  each 
constituent  of  the  microstructure  under  consideration.  This  approach  is  unique  in  that  an 
anisotropic  failure  theory  is  used  on  an  isotropic  matrix  material.  This  complexity  is 
necessitated  by  the  fact  that  the  matrix  failure  behavior  will  be  anisotropic  due  to 
microstructural  geometry.  Consider  a  transversely  isotropic  unidirectional  composite.  If 
all  fibers  were  removed,  leaving .  their  holes,  only  a  matrix  of  Swiss  cheese  would 
remain.  Because  of  the  remaining  microstructure,  macroscopic  failure  of  the  material 
will  be  fundamentally  different  in  axial  versus  transverse  directions.  This  results  in  a 
transversely  isotropic  failure  envelope. 


Anisotropic  failure  behavior  of  in  situ  matrix  material  precludes  using  a  simple  isotropic 
failure  criterion  such  as  maximum  principal  stress.  Also,  recognizing  that  constituents 
typically  have  different  ultimate  strengths  in  tension  and  compression,  each  constituent 
failure  criterion  has  a  tensile  and  compressive  subform. 

As  a  first  approximation,  we  would  like  to  simplify  equation  (4.5) .  Pipes  and  Cole  [24] 
demonstrated  some  of  the  difficulties  in  experimentally  determining  stress  interaction 
terms  such  as  FI2.  Further,  Narayanaswami  [25]  demonstrated  numerically  that  setting 
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the  stress  interaction  term  F/2  to  zero  in  the  Tsai-Wu  quadratic  failure  criterion 
(analogous  to  C12  in  equation  (4.5))  in  plane  stress  analyses,  resulted  in  less  than  10% 
error  for  all  the  load  cases  and  materials  considered.  Hence,  we  set  Cn  equal  to  zero. 
Tsai  and  Wu  identify  the  linear  terms  in  equation  (4.7)  as  necessary  to  account  for 
internal  stresses.  Internal  stresses  are  accounted  for  in  the  formulation  of 
Multicontinuum  Theory  so  the  linear  terms  are  eliminated  from  equation  (4.5).  If  later 
analytical  comparisons  against  experimental  results  do  not  provide  a  satisfactory 
correlation,  these  terms,  along  with  the  term  Cn,  will  be  reexamined  for  their  potential 
contributions. 

Noting  the  above,  the  general  form  for  a  stress  interactive  failure  criterion,  after  changing 
to  a  consistent  coefficient  notation,  is  given  by 

AXI\  +  AJl  +  A3I3  +  A4I4  =  1 .  (4.8) 


4.2.1.  Fiber  Failure 

Most  fibers  used  for  composite  reinforcement  have  greater  transverse  strengths  than  the 
matrices  commonly  used  with  them.  Hence,  we  assume  that  transverse  failure  of  these 
composites  is  matrix  dominated.  Based  on  this  assumption,  we  set  A 2  and  A3  equal  to 
zero  in  equation  (4.8)  as  their  associated  stress  invariants  involve  transverse  normal 
stresses.  The  fiber  failure  criterion  reduces  to: 

A{(l{)2  +  A{  l {  =1.  (4.9) 
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To  determine  coefficients  for  each  stress  term,  we  solve  equation  (4.9)  considering 
individual  load  cases  of  pure  in-plane  shear,  tension,  and  compression.  For  the  case  of 
in-plane  shear  load  only  (cri/=  0),  we  find 

Af  =  1 . . 

4  ter 

For  the  case  of  tensile  load  only  (gi/>  0;  <Jn^  ~  0), we 


For  the  case  of  compression  load  only  (cr i/<  0;  on*  -  0),  we  find 

-Af  - _ I _ 

1  ~ur 

The  criterion  for  fiber  failure  can  now  be  expressed  as: 

±A{  {l{  )2  +  A{l{  =  1 .  (4.10) 

The  ±  symbol  indicates  that  the  appropriate  tensile  or  compressive  ultimate  strength  value 
is  used  depending  on  the  constituent’s  stress  state.  All  of  the  fiber  failure  parameters  can 
be  experimentally  determined  using  simple  uniaxial  or  shear  tests. 
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4.2.2.  Matrix  Failure 

To  determine  the  coefficients  of  equation  (4.8)  for  matrix  failure,  we  first  solve  the 
equation  considering  load  cases  of  pure  in-plane  and  transverse  shear.  For  the  case  of 
transverse  shear  only  (07 /  m=022  m~<?33  m-^i2m=  0),  we  find 


For  the  case  of  in-plane  shear  only  (07  /  m=(722  m=<J33  m=023 m-  0),  we  find 


A 


m 

4 


Since  most  fibers  used  for  composite  reinforcement  have  greater  longitudinal  strengths 
than  the  matrices  commonly  used  in  conjunction  with  them,  we  assume  that  the 
longitudinal  failure  of  these  composites  is  fiber  dominated.  Based  on  this  assumption 
(and  some  numerical  sensitivity  studies),  we  set  Ai  equal  to  zero.  Incorporating  these 
results  into  (4.8)  gives 


/3m  + 


(4.11) 


To  determine  A™,  we  consider  the  case  of  transverse  tensile  load  only  (023  m=<Ji2  m=0; 
(cT22m+&33m)>0)  and  find 
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Note  that,  while  a  pure  transverse  (one-dimensional)  load  on  a  composite  lamina  results 
in 

1  :  ^33  =  ® » 


the  constituents  experience  a  fully  three-dimensional  stress  state.  Likewise,  for  the  case 
of  a  pure  transverse  compressive  load  (oy?  m=<Ji2  m=0;  (<722  m+<733  )  <  0) 


The  criterion  for  matrix  failure  can  now  be  expressed  as 


±A?  (/”  )2  +  A3mI3m  +  A”  I”  =  1 .  (4.12) 

All  of  the  matrix  failure  parameters  can  be  experimentally  determined  using  simple 
uniaxial  or  shear  tests. 

4.3.  Material  Point  Failure 

MCT’s  fundamental  concept  of  average  field  values  at  a  continuum  point  is  also  applied 
to  material  properties.  There  is  no  variation  in  material  properties  across  a  constituent  at 
a  continuum  point.  Therefore,  when  a  constituent  material  fails,  it  fails  absolutely  within 
the  point,  without  directionality,  and  loses  the  ability  to  sustain  any  load.  Furthermore,  at 
every  point  in  the  composite,  both  fiber  and  matrix  are  assumed  to  exist.  Therefore, 
when  one  constituent  fails  the  other  constituent  remains. 
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As  discussed  previously,  an  interesting  ramification  of  this  type  of  constituent  based 
failure  is  that  even  though  the  matrix  material  is  isotropic  (Table  4.1)  in  a  neat  or  bulk 
resin  form,  its  in  situ  failure  properties  will  be  anisotropic.  Take  for  example  the 
E-glass/vinylester  composite  used  in  this  research.  From  an  MCT  perspective,  a 
composite  with  failed  fibers  is  composed  of  only  matrix.  Yet  the  “composite”,  composed 
of  matrix  only  results  in  transversely  isotropic  properties  as  determined  by 
micromechanics  (Table  4.2). 

Table  4.1  Vinylester  matrix  properties. 


E  (GPa) 

V 

4.66 

0.292 

Table  4.2  E-glass/vinylester  composite  properties  with  fiber  failure. 


Eh  (GPa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa) 

V12 

FVF 

2.32 

1.11 

0.610 

0.403 

0.290 

51% 
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SECTION  5. 

DEVELOPING  MATERIAL  PROPERTIES 

5.1.  Experimental  Determination  of  Composite  Material  Properties 
The  MCT  77  finite  element  program  is  capable  of  modeling  nonlinear-elastic  composite 
laminates  formed  by  stacking  lamina  composed  of  unidirectional  reinforcing  fibers  in  a 
polymer  matrix.  Theory  and  software  assume  transversely  isotropic  lamina. 

The  key  philosophy  in  conducting  a  MCT  77  failure  analysis  is  to  “tune  the  finite 
element  code  to  materials  composing  the  structure.  Definition  of  a  “new  material 
requiring  its  own  characterization  must  include  not  only  a  composite  with  different 
constituents  but  also  composites  formed  with  different  fiber  volume  fractions  and 
manufacturing  processes.  Stress-strain  curves,  especially  those  with  significant  non¬ 
linear  behavior,  the  moduli  derived  from  them,  and  ultimate  failure  strengths  are  all 
important  parameters  in  characterizing  material  behavior. 

A  suite  of  six  tests  conducted  for  the  most  part  on  unidirectional  laminates  is  necessary  to 
fully  characterize  a  material  for  failure  analysis.  These  tests  are  used  to. 

1)  Determine  nonlinear  stress-strain  parameters. 

2)  Determine  five  independent  composite  elastic  constants:  En,  E22,  G]2,  G23,  and  v12. 

3)  Determine  six  ultimate  composite  strengths: +Sn,  S22,'Sn,  S22,  Si2,  and  S23. 

The  suite  of  material  tests  used  in  this  research  and  their  relevant  ASTM  standards  are 
listed  in  Table  5.1.  It  should  be  noted  that  ASTM  test  methods  recommend  a  minimum 
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of  five  specimens  be  tested  for  statistical  significance  [26,27].  MCT_77  software  does 
not  differentiate  between  tensile  and  compressive  elastic  constants.  Arbitrarily,  tensile 
elastic  constants  were  used  in  all  MCT_77  analyses  contained  within  this  research. 


Table  5.1  ASTM  test  suite  for  characterizing  composite  (lamina)  materials. 


ASTM  STANDARD 

TYPE  OF  TEST 

CONSTANTS  DETERMINED 

D3039 

Longitudinal  Tension 

+  T~I  +  +  Q 

£\U  V12j  ^11 

tt 

Transverse  Tension 

+-j~i  + 

£>22,  V2i,  022 

D3410 

Longitudinal  Compression 

Ell,  v12>  S„ 

it 

Transverse  Compression 

E22?  V2l,  S22 

D3518 

Longitudinal  Shear 

G12,  S12 

D5379 

Transverse  Shear 

G23,  S23 

5.2.  Determining  Constituent  Material  Properties 

MCT_77’s  ability  to  calculate  accurate  constituent  stress  and  strain  fields  is  dependent  on 
constituent  elastic  constants  derived  from  experimentally  determined  composite  values. 
Further,  MCT_77’s  ability  to  execute  realistic  failure  analysis  is  dependent  on  accurate 
values  for  constituent  strengths  also  derived  from  experimentally  determined  composite 
values.  The  relationship  between  composite  (macro)  and  constituent  (micro)  elastic 
constants  is  developed  using  a  finite  element  idealization  of  a  fiber  imbedded  in  a  matrix 
known  as  a  micromechanics  model.  The  finite  element  micromechanics  model  used  in 
this  research  was  advanced  by  Gamich  [13]  which  contains  discussion  of  its 
development.  Only  major  components  of  the  model  will  be  summarized  here. 

The  micromechanics  model  is  based  on  an  assumption  of  uniform  hexagonal  fiber 
packing  within  the  lamina’s  matrix  (Fig.  5.1).  A  unit  cell,  representative  of  the  repeating 
microstructure,  is  extracted  from  a  region  bounded  by  symmetry  lines.  Unit  cell 
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geometry,  fiber  volume  fraction,  and  boundary  conditions  are  used  to  define  the  finite 
element  model  (Fig.  5.2). 


The  unit  cell  is  based  on  a  generalized  plane  strain  assumption  in  the  fiber  direction  but  is 
fully  three-dimensional.  The  cell  is  modeled  with  ANSYS  finite  element  software  using 
an  ANSYS®  scripting  language  that  allows  material  properties  and  fiber  volume  fraction 
to  be  varied  as  required.  The  model  uses  528  eight-node  brick  (SOLID45)  elements  for  a 
total  of  3444  degrees  of  freedom.  Boundary  conditions  necessary  to  enforce 
compatibility  of  unit  cell  boundaries  with  adjacent  unit  cells  are  generated  automatically. 
Four  linear  elastic  load  cases  are  solved  (longitudinal  tension,  transverse  tension, 
transverse  shear  and  longitudinal  shear)  to  determine  and  verify  five  independent  elastic 
constants  (En,  E22,  G12,  G23,  V12)  for  transversely  isotropic  composite  lamina. 
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Fig.  5.2  Finite  element  model  of  a  unit  cell. 


Ideally,  development  of  a  consistent  set  of  composite  elastic  constants  would  begin  by 
acquiring  elastic  constants  for  each  constituent.  Bulk  matrix  material  properties  are 
typically  isotropic  and  could  be  determined  directly  by  standard  (ASTM)  material  testing 
methods.  Reinforcement  properties  are  more  difficult  to  determine  directly  using 
standard  testing  techniques  due  to  their  small  diameters  and  fibrous  nature.  Generally, 
fiber  properties  are  inferred  or  “backed  out”  from  results  of  composite  lamina  tests.  In 
situ  material  constants  for  commonly  used  reinforcing  fibers  and  polymer  matrices  are 
available  in  published  literature  [28]. 

As  reported  elsewhere  [7,5],  the  in  situ  constituent  material  constants  are  generally 
different  than  those  for  the  bulk  materials.  These  differences  can  be  illustrated  by 
examining  matrix  in  situ  and  bulk  shear  behavior.  For  example,  a  composite  lamina’s 
response  to  a  shear  load,  whether  in-  or  out-of-plane,  is  dominated  by  load  transfer  within 
the  matrix.  The  experimentally  determined  shear  stress-strain  curves  for  bulk  vinylester 
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resin  and  a  E-glass/vinylester  lamina  in  Fig.  5.3  not  only  have  different  magnitudes  but 
also  different  shapes.  The  neat  resin  curve  is  nearly  linear  where  composite  curves  show 
significant  nonlinearity. 


Fig.  5.3  Comparison  of  composite/matrix  shear  behavior. 

As  a  result  of  differing  bulk  versus  in  situ  behavior,  all  constituent  material  constants  for 
MCT_77  analyses  were  backed  out  from  experimentally  determined  composite  values  via 
the  micromechanics  model.  Bulk  constituent  material  properties  provided  good  starting 
values  for  the  micromechanics  model  in  an  iterative  approach  to  determine  in  situ 
properties.  MCTJ77,  in  conjunction  with  the  ANSYS®  pre-processor,  can  run  the 
micromechanics  model  in  reverse  so  that  constituent  properties  consistent  with 
experimentally  determined  composite  properties  can  be  determined  in  an  automated 
fashion. 
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5.3.  Constitutive  Model  for  Changing  Shear  Properties 

A  majority  of  composite  materials  in  use  today  have  organic  matrices  that  produce 
significant  nonlinear  shear  stress-strain  behavior  (see  for  example  Fig.  5.3).  Failure  of  a 
structure  can  be  defined  in  terms  of  strain,  displacement,  or  stress.  Hence,  to  completely 
capture  the  response  of  a  structure  to  a  generally  applied  load,  both  the  stress  and  strain 
states  need  to  be  determined.  Failure  to  model  nonlinear  shear  behavior  can  lead  to  poor 
approximations  of  a  composite’s  response  to  load  as  shown  Fig.  5.4.  Nonlinear  behavior 
can  also  be  found  in  transversely  loaded  lamina.  However,  transverse  nonlinearity  occurs 
over  a  much  smaller  strain  range  due  to  low  ultimate  strength  in  that  direction  and  thus 
can  be  neglected  in  most  structural  analyses. 

As  previously  mentioned,  the  success  of  an  MCT-based  analysis  depends  in  large  part  on 
a  consistent  relationship  between  composite  and  constituent  material  properties. 
Analyses  having  significant  shear  require  incremental  application  of  the  load  and  a 
corresponding  application  of  an  appropriate  tangent  shear  modulus  to  account  for 
nonlinear  structural  response  over  the  entire  load  range.  Micromechanics  arguments 
suggest  that  maintaining  consistent  composite-constituent  material  relationships  when 
nonlinear  composite  shear  moduli  exists  also  requires  the  constituents’  shear  moduli  to 
change.  Changes  to  constituent  properties,  in  turn,  affect  other  composite  elastic 
constants  besides  shear  modulus  since  they  are  linked  through  common  constituents. 
Assuming  shear  behavior  of  the  fiber  remains  linear  elastic,  the  matrix  shear  modulus 
must  be  adjusted  to  maintain  consistency  with  the  composite  elastic  constants. 
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Fig.  5.4  MCT_77  laminate  analysis  of  boron/5505  [±30Jswith  and  without  nonlinear 
shear  effects. 


Researchers  have  applied  many  different  constitutive  models  to  address  changing 
material  constants  under  shear  and  combined  load  states.  Most  appear  to  rely  on 
plasticity  [7,29,30,31]  or  nonlinear-elasticity  [32,33,34]  theories.  While  these 
constitutive  models  have  resulted  in  good  approximations  to  composite  nonlinear  shear 
behavior,  others  [13,35,36]  have  shown  that  nonlinear  polymeric  behavior  is  viscoelastic 
in  nature. 


Figures  5.5-5.7  present  hypothetical  stress-strain  curves  for  elastic-plastic,  viscoelastic, 


and  nonlinear-elastic  behaviors.  A  linear  region  terminated  by  a  transition,  or  yield  point, 
to  nonlinear  behavior  characterizes  the  elastic-plastic  curve,  shown  in  Fig.  5.5. 
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Fig.  5.5  Elastic-plastic  stress-strain  curve  for  a  hypothetical  material. 


Upon  unloading,  the  elastic-plastic  curve  reverts  to  its  previous  linear  behavior  but  does 
not  return  to  zero  strain.  The  residual  strain  is  the  result  of  inelastic  plastic  flow. 


The  nonlinear-elastic  curve,  shown  in  Fig.  5.6,  is  nonlinear  throughout  its  load  history 
and  unloads  back  down  the  loading  path  to  zero  strain  if  no  damage  in  the  material  has 
occurred.  The  unloading  path  shown  in  Fig.  5.6  is  offset  slightly  from  the  loading  path 
for  clarity. 
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5.6  Nonlinear-elastic  stress-strain  curve  for  a  hypothetical  material 


Fig.  5.7  Viscoelastic  stress-strain  curve  for  a  hypothetical  material. 
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Viscoelastic  behavior,  shown  in  Fig.  5.7,  is  time  dependent,  whereas  elastic-plastic  and 
nonlinear-elastic  behaviors  are  not.  The  stress-strain  curve  is  nonlinear  throughout  its 
load  history  and  returns  to  zero  strain  at  zero  load  if  no  damage  in  the  material  has 
occurred. 

Figures  5.8  and  5.9  show  experimental  stress-strain  curves  for  two  different  composite 
materials.  The  shape  of  these  curves  support  the  supposition  that  composite  material 
behavior  is  viscoelastic.  However,  from  Figures  5. 5-5.7  it  is  apparent  that  the  three 
constitutive  models,  elastic-plastic,  viscoelastic,  and  nonlinear-elastic,  will  produce 
similar  stress-strain  curves  for  monotonic  loading.  Hence,  for  the  purpose  of  a  failure 
analysis  using  stress-based  failure  criteria,  each  of  the  constituent  models  would  produce 
essentially  identical  failure  results  under  quasi-static  loading  conditions.  Application  of 
Multicontinuum  Theory  using  viscoelasticity  (or  plasticity)  constitutive  models  would  be 
far  more  complex  than  using  a  nonlinear-elastic  model.  Further,  for  the  research 
considered  herein,  unloading  of  the  composite  was  not  a  consideration.  Hence,  a 
nonlinear-elastic  constitutive  model  relating  changes  in  elastic  constants  due  to  changing 
composite  shear  modulus  was  developed.  Linear  elasticity  was  used  to  model  pure 
tension  and  compression  behavior. 
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MCT_77’s  nonlinear-elastic  constitutive  model  was  predicated  on  the  following: 

1 .  Matrix  compressibility  does  not  change  with  stress  state  (constant  bulk  modulus). 

2.  Fiber  elastic  constants  do  not  change  with  stress  state. 

3.  To  correctly  respond  to  the  known  and  limiting  cases  of  pure  in-plane  and  pure  out- 
of-plane  shear  loads,  shear  moduli  derived  from  curve  fits  of  experimentally 
determined  composite  shear  stress-strain  curves  were  used  to  define  nonlinear 
behavior. 

4.  In  cases  where  two  in-plane  or  in-plane  and  out-of-plane  shear  loads  are  present,  the 
shear  contributing  the  largest  amount  of  deviatoric  strain  energy  determines  the 
value  of  the  composite  tangent  shear  modulus  used  in  the  constitutive  model. 


MCT  77  uses  a  three-term  exponential  series  to  fit  both  in-  and  out-of-plane 
experimental  shear  stress-strain  curves.  Curve  fits  of  shear  stress-strain  data  have  the 
form 


t  =  C0  +  C,e 


+  C2e 


(«2  r) 


(5.1) 


where  Q  and  a,  are  curve  fit  parameters  and  y  is  engineering  shear  strain.  Nonlinear 
regression  is  used  to  fit  the  five  equation  parameters  to  data  generated  by  experimental 
tests.  Strain  dependent,  tangent  shear  moduli  are  computed  from  the  first  derivative  of 
equation  (5.1)  for  use  during  a  finite  element  analysis.  Appendix  D  provides  a  detailed 
discussion  of  the  nonlinear  regression  procedure  used. 


During  an  MCT_77  analysis,  recalculation  of  all  composite  and  matrix  material  constants 
is  initiated  by  changes  in  composite  shear  strain  values.  The  relationships  between 
composite  and  matrix  material  constants  are  developed  in  an  unconnected  pre-processing 
procedure  by  varying  matrix  properties  in  a  prescribed  manner  and  using  micromechanics 
to  derive  new  composite  properties.  The  procedure  varies  the  matrix  shear  modulus  from 
0%-100%  (usually  in  increments  of  10%)  of  its  initial  in  situ  value  as  determined  from 
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material  tests.  An  assumption  of  a  constant  bulk  modulus  provides  the  second 
independent  elastic  constant.  Discrete  increments  in  matrix  elastic  constants  and  constant 
fiber  elastic  constants  form  input  to  the  micromechanics  model  for  calculating 
corresponding  composite  constants.  A  curve  for  each  elastic  constant,  normalized  with 
respect  to  the  composite  tangent  shear  modulus,  was  fitted  through  the  values  determined 
at  each  data  point.  Quadratic  equations  for  four  curves,  E22,  V12,  G23,  and  Gm,  referenced 
to  changes  in  composite  G12  and  four  curves,  E22,  V12,  G12,  and  Gm,  referenced  to  changes 
in  composite  G23  were  incorporated  into  the  MCT  77  program.  An  example  of  a 
quadratic  shear  relationship  is  shown  in  Fig.  5.10. 


All  quadratic  curve  fits  of  the  micromechanical  data  were  excellent,  resulting  in  R  values 
of  0.97  or  better.  Appendix  A  contains  the  complete  set  of  quadratic  curve  fits  for 
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boron/epoxy.  Appendix  B  for  E-glass/vinylester,  and  Appendix  C  for  carbon/epoxy 
composites  used  in  this  research. 

5.4.  Determining  Constituent  Ultimate  Strengths 

The  MCT  failure  criterion  developed  previously  required  determining  five  failure 
parameters,  A2f,  A4{  A2  m,  Aj  and  A4  m.  These  parameters  are  functions  of  ultimate 
constituent  strengths,  ±Suf,  Sl2f,  ±S22  m,  "Sjj  m,  S22  and  S22  m.  Constituent  strengths  were 
derived  from  experimentally  determined  composite  lamina  ultimate  strengths. 
Determining  which  constituent  precipitates  composite  failure  is  necessary  for  establishing 
accurate  constituent  failure  values.  Identifying  the  constituent  that  precipitates  failure  in 
longitudinal  and  transverse  lamina  tension  and  compression  tests  is  intuitive  and 
straightforward,  i.e.,  fiber  failure  for  longitudinal  loads  and  matrix  failure  for  transverse 
loads.  Identifying  the  constituent  leading  to  shear  failure  is  more  problematic,  as  non- 
catastrophic  matrix  and  fiber  damage  begins  well  before  ultimate  composite  strength  is 
achieved  [37].  Good  initial  values  for  constituent  shear  strengths  could  be  obtained  from 
an  MCT  analysis  of  experimentally  determined  ultimate  composite  shear  strengths. 

In  an  effort  to  more  rigorously  determine  the  ultimate  constituent  shear  strengths,  a 
procedure  using  nonlinear  regression  analysis  of  load  cases  involving  shear  was 
developed.  Data  from  off-angle,  balanced,  symmetric  laminates,  [±0]s,  provided  an 
excellent  basis  for  determining  optimized  constituent  failure  parameters.  These  laminates 
produced  varying  degrees  of  combined  shear  and  normal  stresses  and  tended  to  fail  in 
modes  that  allowed  analytical  identification  of  the  constituent  precipitating  laminate 
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failure.  Off-angle,  unbalanced,  laminates,  [+6]n  or  [-0]n  also  produced  varying  shear  and 
normal  stresses  but  had  to  be  used  with  caution  because  of  non-uniform  stress  states 
caused  by  warping  of  the  test  specimen  [24]. 


Using  ultimate  composite  failure  strengths  for  each  off-angle  laminate,  MCT_77 
calculated  the  constituents’  transversely  isotropic  stress  invariants  at  the  experimentally 
determined  point  of  composite  failure.  Next,  the  failure  mode  for  each  laminate  tested 
was  tentatively  identified.  As  a  general  rule,  composite  failure  is  fiber  dominated  at 
lower  fiber  orientation  angles  (approximately  20°  and  less)  and  matrix  dominated  at 
higher  fiber  angles.  Visually  determining  failure  modes  from  test  specimens  is  difficult, 
so  it  becomes  more  convenient  to  do  analytically.  As  mentioned  previously,  good  initial 
guesses  for  constituent  ultimate  shear  strengths  using  MCT_77  can  be  derived  from 
experimentally  determined  unidirectional  lamina  ultimate  shear  strengths  or  taken  from 
the  available  literature.  The  data  were  sorted  by  failure  mode  (fiber  or  matrix).  Using  the 
failure  criterion  and  nonlinear  regression  analysis,  a  best  fit  was  determined  for  the 
failure  parameter  Aj  m  and  A4  m  for  the  matrix  and  A 4 ^  for  the  fiber.  After  the  first  trial 
optimization,  updated  values  for  As  m,  A4  m,  and  A4  f  were  used  to  recheck  the  failure 
modes.  If  modes  changed,  another  iteration  was  required.  Once  failure  parameters  were 
“tuned”  to  off-angle  laminate  results,  MCT_77  failure  analyses  on  more  complicated 
composite  structural  laminates  were  conducted.  Appendix  E  provides  a  detailed 
discussion  of  the  optimization  procedure  used. 
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SECTION  6. 

MCT  77  FINITE  ELEMENT  CODE 


6.1.  General  Capabilities  of  MCT_77 

The  goal  of  the  research  contained  within  this  dissertation  was  to  develop  a  finite- 
element-based  tool  for  the  failure  analysis  of  practical  composite  structural  laminates.  To 
that  end,  a  “building  block”  philosophy  was  used  to  develop  the  finite  element  software: 


1 .  Discretize  (“mesh”)  the  structure  using  standard  FEA  techniques  and  commonly 
available  graphical  software. 

2.  Exercise  judicious  element  choices.  Three-dimensional  solid  elements  are  preferred 
in  most  cases  because  they  produce  a  complete  stress  tensor.  The  three- 
dimensional  capability  penalizes  the  designer  by  requiring  re-meshing  of  the 
structure  with  any  changes  in  laminate  thickness.  Shell  elements  are  overwhelming 
favored  by  designers  because  they  allow  element  thickness  to  be  specified  without 
changing  the  original  structural  element  mesh.  Shell  elements  provide  only  five  out 
of  the  six  stress  tensor  components  because  they  are  formulated  on  a  modified 
plane-stress  assumption.  The  normal  through-thickness  stress,  important  in 
delamination,  is  assumed  to  be  zero. 

3.  Use  a  layering  capability  in  each  element  to  breakdown  the  composite  laminate  into 
individual  lamina. 

4.  Breakdown  the  lamina  into  constituents  (fiber  and  matrix)  using  MCT. 

5.  Conduct  failure  analysis  at  the  constituent  level  where  structural  damage  initiates. 


MCT_77  is  a  conventional  displacement-based,  finite  element  code.  What  sets  MCT  77 
apart  from  other  finite  element  codes  is  an  implementation  of  Multicontinuum  Theory 
allowing  constituent  stress  and  strain  fields  to  be  extracted  from  those  of  the  composite 
material.  Material  behavior  is  limited  to  nonlinear  elastic.  Inelastic  behavior  due  to 
damage  is  accommodated  through  a  progressive  failure  algorithm.  MCT_77  was  written 
in  FORTRAN  77  for  maximum  portability  between  computer  platforms.  The  code  was 
extensively  commented  and  written  in  a  modular  fashion  for  readability  and  to  facilitate 
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upgrades.  Appendix  F  provides  a  detailed  explanation  of  the  MCT_77  program  and  its 
inputs. 


6.1.1.  Element  Library 

MCTJ77  currently  has  a  four-element  library  consisting  of  the  following: 


1 .  QUAD4:  Four-node,  two-dimensional,  2  degree-of-ffeedom  (DOF),  plane 
stress/strain,  bilinear,  isoparametric,  quadrilateral  element. 

2.  QUAD8:  Eight-node,  two-dimensional,  2  DOF,  plane  stress/strain,  quadratic 
(serendipity)  isoparametric,  quadrilateral  element. 

3.  BRICK8:  Eight-node,  three-dimensional,  3  DOF,  trilinear,  isoparametric 
hexahedral  element. 

4.  SHELL8:  Eight-node,  two-dimensional,  5  DOF,  modified  plane  stress,  quadratic 
(serendipity),  isoparametric  shell  element. 


Currently,  elements  cannot  be  mixed  within  a  MCT  77  analysis  nor  can  they  be 
degenerated  to  triangular  or  wedge  shapes. 


The  linear  elements,  QUAD4  and  BRICK8,  have  the  option  of  “incompatible”  (or 
“nonconforming”)  element  formulations  [38].  This  formulation  uses  shape  functions 
incorporating  internal  degrees-of-ffeedom  permitting  displacements  to  be  discontinuous 
across  adjacent  element  boundaries,  hence,  "incompatible".  The  formulation  allows 
element  boundaries  to  assume  curved  shapes  in  a  manner  similar  to  higher  order 
quadratic  elements  but  without  the  additional  nodes  and  associated  increase  in 
computational  and  storage  requirements.  In  the  presence  of  shear,  incompatible  element 
formulation  results  in  more  stable  element  behavior  and  greater  accuracy  than  regular 
linear  elements. 
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All  MCT  77  elements  have  the  ability  to  model  lamina  layering.  The  user  specifies  total 
laminate  thickness  and  number  of  lamina.  Currently,  MCT_77  assumes  a  constant 
lamina  thickness.  Formulation  of  lamina  layer  stiffness  is  the  same  as  the  element  that 
contains  it.  SHELL8  element  layers  are  under  a  modified  plane  stress  assumption 
(ojj=0)  and  use  2x2x2  (in-plane  element  x  direction  by  in-plane  element  y  direction  by 
out-of-plane  element  z  direction)  Gauss  integration  within  each  layer.  QUAD4  and 
QUAD8  layers  are  under  full  plane  stress  assumption  (cr?j= 073=0^=0)  and  use  2x2x1 
Gauss  integration  within  each  layer.  BRICK8  layers  produce  a  full  stress  tensor  and  use 
2x2x2  Gauss  integration  within  each  layer. 

6. 1 .2.  Boundary  Conditions 

MCT_77  is  specifically  designed  to  support  boundary  conditions  that  typically  arise  in 
structural  analysis.  MCT_77  supports  any  admissible  combination  of  the  following 
prescribed  loads: 

a)  Nodal  forces. 

b)  Uniform  temperature. 

c)  Nodal  temperatures. 

d)  Element  pressures. 

e)  Body  forces  arising  from  translational  or  rotational  accelerations. 

MCT_77  supports  any  admissible  combination  of  the  following  prescribed 
displacements: 

a)  Homogeneous  nodal  displacements. 

b)  Non-homogeneous  nodal  displacements. 

c)  Nodal  constraint  equations. 
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MCT  77  can  apply  loads  and  non-homogeneous  nodal  displacements  as  a  single  step  or 
arbitrary  ramp  function.  The  incremental  nature  of  ramping  loads  from  zero  to  a  target 
value  lends  itself  to  the  nonlinear  behavior  characteristic  of  structural  failure  analysis. 

Constraint  equations  implemented  in  MCT_77  allow  the  user  to  define  complex 
relationships  between  degrees-of-freedom  of  different  nodes.  For  example,  constraint 
equations  are  useful  in  applying  boundary  conditions  representing  planes  of  symmetry. 
Constraint  equations  have  the  form: 

C,  *  DOF?  +C2  *  DOF]  +...  +  C„  *  DOF”  =  C 

where  C,  are  constants  and  DOF"  are  degree  of  freedom  i  of  node  n.  A  simple  constraint 
equation  coupling  the  x  degrees-of-freedom  of  nodes  1  and  2  so  they  move  together 
(equal  magnitude)  would  look  like 

(l)  *  DOF?  +  (- 1)  *  DOF2  =  0 . 

6.2.  ANSYS®  Interfacing 

A  key  component  to  a  successful  finite  element  analysis  is  the  assimilation  and 
interpretation  of  results  by  the  user.  Even  a  modest  size  analysis  can  result  in  a 
voluminous  amount  of  data  describing  the  distribution  and  intensity  of  field  variables 
over  the  structural  domain.  Further,  inputting  a  structure’s  geometry  and  developing  a 
finite  element  mesh  as  a  sequence  of  formatted  numbers  within  a  text  file  can  be  tedious 
and  error  prone. 
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It  is  beyond  the  scope  of  this,  or  perhaps  any,  structural  research  effort  to  attempt 
independent  development  of  a  visual  simulation  program  when  so  many  are 
commercially  available.  A  more  efficient  approach  is  to  write  computer  programs  that 
function  as  an  interface  between  finite  element  input  and  output  information  with  that 
required  by  a  commercially  available  visual  simulation  program. 

To  facilitate  construction  of  finite  element  models  for  MCT_77  analyses,  a  translator, 
ANS2MCT,  has  been  written  allowing  the  designer  to  use  pre-processing  capabilities  of 
the  ANSYS  [39]  finite  element  program  to  create  the  geometry  and  finite  element  mesh 
of  the  structure  being  analyzed.  ANS2MCT  converts  ANSYS®  database  information  to 
MCT_77  readable  format.  If  requested,  the  MCT  77  code  writes  results  files  that  can  be 
read  by  ANSYS®  for  post-processing. 

As  shown  in  Figures  6. 1-6.6,  visualization  can  be  especially  useful  in  progressive  failure 
analyses  when  specific  color-coding  designates  damaged  regions  within  a  laminate.  This 
allows  the  user  to  step  through  the  load  history  and  watch  material  failure  initiate  and 
propagate  throughout  the  laminate.  A  user/designer  can  then  modify  the  material 
architecture  of  the  structure  in  a  more  precise  manner  tending  to  optimize  the  design 
goals  of  maximum  structural  strength  or  stiffness  at  minimum  weight  and  cost.  The  color 
codes  in  Figures  62-6.6  use  blue  to  indicate  undamaged  composite  material,  green  to 
indicate  composite  damage  due  to  matrix  failure,  and  red  to  indicate  composite  damage 
due  to  fiber  failure.  Section  6.5  discusses  an  MCT  based  progressive  failure  analysis 
approach  in  detail. 
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Fig.  6.1  FE  idealization  of  a  cruciform  specimen  under  biaxial  tension-tension  load. 


Fig.  6.2  Composite  damage  in  the  cruciform  specimen  at  load  step  9. 
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Fig.  6.3  Composite  damage  in  the  cruciform  specimen  at  load  step  10. 


Fig.  6.4  Composite  damage  in  the  cruciform  specimen  at  load  step  22. 
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Fig.  6.5  Composite  damage  in  the  cruciform  specimen  at  load  step  23. 


ME  BOTH 


Fig.  6.6  Composite  damage  in  the  cruciform  specimen  at  load  step  24. 
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6.3.  Incorporation  of  Multicontinuum  Theory 

MCT  equations,  developed  in  Section  3,  are  implemented  within  the  MCT_77  code  in 
two  subroutines.  The  first  subroutine  calculates  the  [A]  and  {a}  (see  equation  3.13) 
matrices  for  all  composite  materials  defined  by  the  user.  The  second  subroutine  extracts 
constituent  strains  from  composite  values  and  computes  constituent  stresses. 
Calculations  within  these  routines  follow  their  derivation  and  are  applied  in  a 
straightforward  manner. 

For  a  linear  elastic  analysis  in  which  material  properties  do  not  change,  an  MCT 
decomposition  of  composite  strains  and  stresses  could  be  accomplished  in  a  single  post¬ 
processing  step.  The  brevity  and  simplistic  character  of  the  MCT  specific  subroutines 
offers  the  potential  for  a  relatively  quick  and  trouble-free  incorporation  as  user-defined 
subroutines  available  as  part  of  most  commercial  finite  element  codes.  However,  for 
nonlinear  progressive  failure  analysis,  complete  integration  of  these  subroutines  with  the 
bulk  FEA  code  is  computationally  advantageous. 


6.4.  Nonlinear  Constitutive  Model  for  Changing  Shear  Properties 
MCTJ77  execution  of  the  constitutive  model  developed  in  Section  5  requires  the  load  to 
be  applied  as  a  step-wise  ramp.  When  shear  stresses  are  present,  MCT  77  computes  a 
composite  tangent  modulus  from  the  appropriate  shear  stress-strain  curve  fit  information. 
Using  the  tangent  shear  modulus,  the  remaining  elastic  constants  for  both  composite  and 
constituents  are  calculated  from  equations  relating  them  to  change  in  shear  modulus  (See 
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Figures  A.6-A.9,  B.7-B.14,  and  C.7-C.14).  Element  stiffness  is  recalculated  based  on 
updated  material  properties.  A  gradual  softening  of  the  structure  causes  a  load 
equilibrium  imbalance  which  is  resolved  by  Newton-Raphson  iterations  as  described  in 
section  6.5.  Once  equilibrium  is  achieved,  the  process  is  repeated  for  subsequent  load 
steps  until  the  maximum  load  is  reached. 

Most  of  the  nonlinear  constitutive  model  is  incorporated  as  an  algorithm  within  a  single 
MCT_77  subroutine.  When  an  analysis  incorporating  nonlinear  shear  is  requested,  the 
algorithm  is  executed  at  every  Gauss  point  within  the  structure.  Executing  a  small 
collection  of  algebraic  equations  to  produce  macro-  and  micromechanical  material  data  is 
several  orders  of  magnitude  faster  than  executing  the  micromechanics  finite  element 
model  itself  and  produces  virtually  identical  results. 


6. 5.  Progressive  Failure  Analysis 

In  the  finite  element  method,  numerical  integration  samples  stress,  strain,  and  material 
values  at  Gauss  quadrature  or  “material”  points.  There  is  a  hierarchy  to  the  material 
numbering  in  the  MCT_77  finite  element  code  that  corresponds  to  severity  of  composite 
damage  at  a  material  point.  Three  composite  material  conditions  or  states,  listed  in 
increasing  severity,  are  defined  as: 

1 .  Undamaged  composite, 

2.  Composite  damaged  by  matrix  failure,  and 

3.  Composite  damaged  by  fiber  failure. 

There  are  four  constituent  states  corresponding  to  the  three  composite  damage  states: 
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1 .  Intact  fiber  (composite  states  1  and  2), 

2.  Failed  fiber  (composite  state  3), 

3.  Intact  matrix  (composite  states  1  and  3),  and 

4.  Failed  matrix,  (composite  state  2  and  3). 

When  either  constituent  fails,  its  moduli  is  reduced  to  a  near  zero  value.  Near  zero  values 
are  used  rather  than  zero  to  avoid  numerical  difficulties.  Matrix  moduli  are  reduced  to 
1%  of  their  original  value.  Fiber  moduli,  which  are  typically  one  to  two  orders  of 
magnitude  larger  than  matrix  moduli,  are  reduced  by  whatever  percentage  is  required  to 
bring*  damaged  fiber  values  to  the  same  magnitude  as  damaged  matrix  so  that  near  zero 
stiffness  values  are  the  same  for  both  constituents.  In  practice,  a  failed  fiber  constituent 
reduces  the  load  carrying  capacity  of  the  composite  so  severely  that  it  can  be  assumed 
that  both  constituents  have  failed.  Since  all  constituent  properties,  both  intact  and  failed, 
are  known  a  priori,  micromechanics  can  be  used  to  determine  damaged  composite 
properties  due  to  a  failed  constituent  outside  the  MCT_77  program.  Files  defining 
damaged  and  intact  composite  and  constituent  elastic  constants  form  input  to  the 
MCT_77  program  as  seven  different  material  types:  three  composite  and  four  constituent 
damage  states. 

A  MCT_77  failure  analysis  requires  the  load  to  be  incrementally  applied.  The  composite 
material  damage  state  at  every  Gauss  point  is  stored  for  the  entire  analysis.  Initially, 
composite  material  properties  are  set  to  a  state  1  (undamaged)  condition.  At  each  load 
step,  a  damage  algorithm,  using  the  failure  criteria  formulated  in  Section  4,  checks  each 
Gauss  point  for  constituent  failure  based  on  accumulative  stresses.  If  constituent  failure 
is  detected,  composite  material  properties  are  set  to  either  state  2  or  3  and  all  stresses 
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(composite  and  constituent)  at  that  Gauss  point  are  set  to  zero.  Gauss  point  stresses  are 
recalculated  using  accumulated  strains  and  updated  material  properties. 

Implementation  of  the  Tsai-Wu  failure  criterion  within  MCT_77  is  a  straightforward 
duplication  of  the  MCT  failure  criterion  implementation.  Instead  of  examining 
constituent  stress  states,  the  program  uses  the  composite  stress  state  at  each  material 
point.  If  Tsai-Wu  indicates  composite  failure,  the  material  properties  at  that  point  are  set 
to  the  most  severe  damage  state  (number  3). 

After  all  Gauss  points  in  the  model  have  been  checked  for  material  failure  and  modified 
(if  necessary),  an  internal  (resisting)  load  vector  is  calculated.  Gradual  softening  of  the 
structure  due  to  Gauss  point  material  failures  will  cause  an  equilibrium  imbalance 
between  the  applied  (external)  and  internal  load  vectors.  An  iterative  procedure  within 
each  load  step  calculates  differences  between  external  and  internal  load  vectors,  called 
the  residual  load  vector,  and  applies  it  to  the  structure  as  a  “virtual”  load.  The  net  effect 
is  to  increase  nodal  displacements.  Hence,  Gauss  point  strains  and  stresses  increase  until 
equilibrium  is  restored  and  the  next  load  step  is  applied. 

The  general  process  outlined  above  is  standard  procedure  for  nonlinear  finite  element 
problems.  Particular  solution  methods  can  be  found  in  the  literature  [40]  but  the  most 
common  methods  are  the  Newton-Raphson  and  Modified  Newton-Raphson.  The 
Newton-Raphson  method  recalculates  the  global  stiffness  matrix  each  time  there  is  a 
change  in  material  properties  of  the  structure.  Formulating  the  global  stiffness  matrix  is 
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the  single  most  time  consuming  task  within  the  finite  element  method.  MCT  77  uses  the 
Modified  Newton-Raphson  method  which  recalculates  the  global  stiffness  matrix  only  at 
the  beginning  of  a  load  step  and  reuses  it  during  subsequent  equilibrium  iterations. 
Reusing  the  global  stiffness  matrix  requires  more  iteration  for  equilibrium  but  results  in 
significant  timesaving. 


Equilibrium  is  achieved  when  a  ratio  of  the  Euclidean  (L2)  norms  of  the  residual-to- 
accumulated  external  load  vectors  is  less  than  a  prescribed  tolerance 


(er2)2 

— -  ■'  \  <  FTOL  i  =  1,  total  DOF .  (6.1) 


MCT  77  allows  user  specified  FTOL  values  but  typically  values  on  the  order  of  10~3 
(10"4  is  default)  provide  sufficient  accuracy. 


Structural  failure  is  defined  as  that  point  in  the  load  history  when  the  structure  can  no 
longer  support  the  accumulated  load  and  deflections  begin  to  grow  without  bound. 
Unbounded  growth  is  detected  during  equilibrium  iterations  by  measuring  the  percentage 
change  in  the  Euclidean  norm  of  subsequent  displacement  vectors  referenced  to  the  norm 
of  the  first  iteration  displacement  vector.  The  measure,  taken  within  a  load  step,  is  given 
by 


(ZAP.^-feAP2)2 

(EAP,j)5 


<  DTOL , 


(6.2) 
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SECTION  7. 

COMPARISON  OF  ANALYSIS  VERSUS  EXPERIMENT 

7.1.  Overview 

In  this  section,  results  from  MCT_77  failure  analyses  are  compared  against  experimental 
data  of  laminates  fabricated  from  three  different  materials  and  tested  under  uniaxial  or 
biaxial  load  conditions.  Stress-strain  curves  are  used  to  qualify  MCT_77’s  ability  to 
simulate  laminate  structural  behavior.  As  mentioned  in  Section  4,  the  Tsai-Wu  criterion 
is  the  most  widely  used  stress-interactive  failure  prediction  methodology.  Hence,  it  was 
considered  in  this  research  as  state-of-the-practice  and  used  in  conjunction  with 
experimental  data  to  benchmark  the  accuracy  of  MCI  failure  predictions. 

7.2.  Uniaxial  Loading 
7.2.1.  Boron/5505  Composite 

Experimental  stress-strain  curves  for  laminate  specimens  fabricated  from  boron/5505 
(boron/epoxy)  composites  were  digitized  from  figures  published  by  Petit  and  Waddoups 
[41].  The  authors  note  that  the  stress-strain  curves  presented  were  based  on  a  single  test; 
hence,  there  was  no  statistical  basis  for  the  results.  MCT_77  results  are  based  on 
constitutive  material  properties  backed  out  via  micromechanics  from  experimental  data  as 
reported  in  Appendix  A. 
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Ten  laminates  of  different  architectures  were  tested  under  uniaxial  tension,  uniaxial 
compression,  or  both.  Stress-strain  graphs  for  all  laminates  tested  by  Petit  and  Waddoups 
are  included  in  Appendix  A.  Only  a  few  of  the  more  interesting  results  are  discussed 

here. 


Fig.  7.1  Boron/5505  [0/90] s  laminate  under  uniaxial  tension. 


Fig.  7.1  presents  the  experimentally  and  analytically  simulated  stress-strain  behavior  of  a 
[0/90]s  laminate.  Analytically,  the  matrix  in  the  laminate’s  90°  plies  experienced  tensile 
failure  at  approximately  500  MPa  but  the  laminate  continued  to  load  without  significant 
loss  in  stiffness  to  approximately  700  MPa.  The  experimental  data  appear  to  support  this 
scenario  with  initially  linear  behavior  transitioning  to  nonlinear  behavior  at 
approximately  400  MPa  due  conceivably  to  matrix  damage.  MCT  and  Tsai-Wu  laminate 
strength  predictions  were  almost  identical  and  within  1  %  of  the  experimental  value. 
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Fig.  7.2  Boron/5505  [±45]s  laminate  under  uniaxial  tension. 

In  Fig.  7.2,  MCT  77  captured  the  highly  non-linear  response  of  a  [±45]s  laminate  and 
predicted  an  ultimate  failure  load  within  1%  of  the  experimental  value.  Tsai-Wu  failure 
prediction  was  within  10%.  The  Tsai-Wu  curve  follows  the  MCT  curve  identically  up  to 
failure  and  is  offset  for  clarity.  MCT  analysis  indicated  laminate  failure  was  precipitated 
by  matrix  tensile  failure. 

Currently,  MCT_77  can  accommodate  only  one  modulus  per  direction,  i.e.;  it  cannot 
differentiate  between  tensile  and  compressive  moduli.  The  boron/5505  lamina 
compressive  Young’s  modulus  (227  GPa)  was  approximately  10%  higher  than  the  tensile 
Young’s  modulus  (207  GPa).  MCT  77  material  behavior  was  based  on  the  tensile 
Young’s  modulus  hence  the  stress-strain  curves  in  analyses  of  all  compressive  tests  were 
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more  compliant  than  those  of  the  experimental  data.  In  practice,  an  average  of  the  tensile 
and  compressive  moduli  values  could  be  used  without  a  significant  decrease  in  analysis 
accuracy. 


Fig.  7.3  Boron/5505  [±20]s  laminate  under  uniaxial  compression. 

The  [±20]s  laminate,  shown  in  Fig.  7.3,  experienced  a  linear  response  throughout  the 
compressive  load.  Analytically,  MCT  predicted  laminate  failure  due  to  matrix  tensile 
rupture.  MCT  predicted  a  laminate  failure  load  within  12%  and  Tsai-Wu  within  17%  of 
the  experimental  result.  Again,  the  Tsai-Wu  curve  is  offset  for  clarity. 
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Fig.  7.4  Boron/5505  [(0)i/±45]s  laminate  under  uniaxial  compression. 


The  [(0)3/±45]s  laminate,  shown  in  Fig.  7.4,  experienced  linear  response  throughout  the 
load.  The  MCT  matrix  failure  criterion  does  not  allow  for  matrix  failure  in  the  fiber 
direction  but  the  analysis  indicated  transverse  tensile  matrix  failure  in  the  0°  plies  at 
approximately  75%  of  the  ultimate  load  due  to  a  Poisson  effect.  Any  indication  of  matrix 
failure  in  the  experimental  data  appeared  near  ultimate  load.  Ultimate  laminate  failure 
occurred  when  fibers  in  the  0°  plies  ruptured.  MCT  ultimate  failure  load  prediction  was 
within  10%  of  the  experimental  value.  The  Tsai-Wu  prediction  was  off  by  26%  of  the 
experimental  value  because  it  coincides  with  transverse  lamina  (matrix)  failure  in  the  0° 
plys.  The  Tsai-Wu  curve  is  offset  from  the  MCT  curve  for  clarity. 
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A  summary  of  MCT  and  Tsai-Wu  failure  predictions  for  the  boron/5505  laminates 
considered  is  presented  in  Table  7.1.  A  measure  of  predictive  accuracy  or  “modeling 
bias”  is  calculated  as  (predicted  failure  value/experimentally  determined  value).  Scatter 
in  the  data  is  measured  as  the  coefficient-of-variation  (COV)  about  a  mean  value  and  is 
calculated  as  ((mean  value/standard  deviation)  x  100). 


Table  7.1  Summary  of  failure  loads  for  boron/5505  laminates. 


Laminate 

Exd 

MCT/Exd 

Tsai-Wu/Exo 

[0/90]s  C 

-1814 

0.92 

0.92 

[0/90]s  T 

699 

1.01 

0.99 

[±20]s  C 

-766 

1.12 

1.17 

[±20]s  T 

927 

1.16 

1.32 

[±30]s  C 

-307 

1.32 

1.32 

[±30]s  T 

664 

0.97 

0.99 

[±45]s  T 

245 

1.00 

0.90 

[±60]s  C 

-397 

0.63 

0.66 

[±60]s  T 

122 

0.93 

0.86 

[(y±45]s  c 

-1833 

1.10 

0.74 

[03/±45]s  T 

941 

0.98 

1.20 

[0/90/±45]s  T 

497 

1.02 

0.82 

[0/±60]s  T 

499 

1.04 

0.97 

[653/20/-70]s  T 

136 

1.47 

1.47 

[903/±45]s  T 

138 

1.13 

1.19 

Avg  = 

1.05 

1.03 

COV  = 

17.9% 

22.6% 

Table  7.1  shows  there  is  a  significant  amount  of  scatter  in  both  failure  predictions  for  the 
boron/5505  uniaxial  load  cases.  This  scatter  prevents  a  statistically  meaningful 
comparison  between  MCT  and  Tsai-Wu  criterion  for  failure  prediction  accuracy.  The 
large  COVs  result  from  both  the  lack  of  an  experimentally  determined  failure  value 
representative  of  the  data  sets  due  to  using  only  one  data  point  per  test  and  inaccuracies  in 
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the  failure  criteria.  It  is  statistically  significant,  however,  that  the  Tsai-Wu  criterion 
produced  greater  scatter  about  the  mean  average  value  than  MCT. 

7.2.2.  E-glass/8084  Composite 

Stress-strain  curves  for  laminate  specimens  fabricated  from  E-glass/8084 
(glass/vinylester)  composites  were  experimentally  determined  by  the  University  of 
Wyoming’s  Composite  Materials  Research  Group  from  laminates  provided  by  Seemann 
Composites  Inc. 

MCT_77  results  were  based  on  constitutive  material  properties  backed  out  via 
micromechanics  from  experimental  data  based  on  at  least  five  samples  per  test.  Material 
properties  used  in  the  analysis  are  discussed  in  Appendix  B.  E-glass/8084  lamina  had 
essentially  the  same  modulus  values  in  tension  and  compression.  Ten  laminates  of 
different  architectures  were  tested  under  uniaxial  tension,  uniaxial  compression  or  both. 
Stress-strain  graphs  for  all  the  E-glass/8084  laminates  tested  are  included  in  Appendix  B. 
A  limited  number  of  results  are  discussed  here. 

Seven  unidirectional  lamina  with  increasing  fiber  angle,  relative  to  the  direction  of  the 
applied  tensile  load,  were  tested  to  develop  data  for  failure  criteria  optimization. 
Agreement  between  MCT  failure  predictions  and  off-angle  tension  test  results,  shown  in 
Fig.  7.5,  was  excellent  and  can  be  attributed  to  self-consistency  of  failure  parameters 
generated  by  the  optimization  procedure  discussed  in  Section  5.  Agreement  between 
Tsai-Wu  predictions  and  test  results  was  also  excellent  and  can  be  attributed  to  the  nature 
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of  the  off-axis  test  in  which  any  failure,  matrix  or  fiber,  results  in  catastrophic  failure  of 
the  composite. 


Fig.  7.5  Off-angle  E-glass/8084  lamina  under  uniaxial  tension. 

Fig.  7.6  shows  the  stress-strain  response  for  a  [0/90]s  laminate.  Analytically,  both  MCT 
and  Tsai-Wu  indicated  tensile  failure  in  the  90°  plies  at  approximately  100  MPa.  The 
laminate  continued  to  load  without  a  significant  loss  in  stiffness  to  approximately  350 
MPa.  As  with  the  boron/5505  [0/90]  laminate,  the  MCT  analysis  indicated  transverse 
tensile  matrix  failure  in  the  0°  plies  near  the  ultimate  load  due  to  a  Poisson  effect.  Both 
MCT  and  Tsai-Wu  predict  catastrophic  laminate  failure  after  tensile  rupture  of  the  0° 
plies.  MCT  under-predicted  laminate  strength  by  3%  and  Tsai-Wu  over-predicted 
laminate  strength  by  1 1%  of  the  experimental  value. 
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The  quasi-isotropic  [0/±45/90]s  laminate,  shown  in  Fig.  7.7,  had  a  tensile  load  response 
reminiscent  of  that  seen  in  the  [0/90]s  laminate  as  matrix  tensile  failure  occurred 
analytically  in  the  90°  plies  at  approximately  80  MPa.  Additional  matrix  tensile  failure  in 
the  ±45°  plies  was  predicted  prior  to  catastrophic  tensile  fiber  failure  in  the  0°  plies. 
Tsai-Wu  predicted  only  tensile  failure  in  the  90°  plies.  MCT  and  Tsai-Wu  stress-strain 
results  were  in  excellent  agreement  with  experimental  data  and  predicted  ultimate  failure 
loads  within  13%  of  the  experimental  value. 


Fig.  7.8  E-glass/8084  [0/90/±45]s  laminate  under  uniaxial  compression. 

The  quasi-isotropic  [0/±45/90]s  laminate,  shown  in  Fig.  7.8,  had  a  linear  compressive 
load  response  to  about  -180  MPa.  At  that  point  both  the  MCT  analysis  and  the  nonlinear 
behavior  of  the  experimental  data  suggested  matrix  compressive  failure  initiated  in  the 
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90°  plies.  The  Tsai-Wu  prediction  was  linear  to  a  failure  load  21%  below  the 
experimental  value.  MCT  predicted  an  ultimate  load  5%  below  the  experimental  value. 


A  summary  of  MCT  and  Tsai-Wu  failure  predictions  for  the  E-glass/8084  laminates 
analyzed  is  presented  in  Table  7.1.  There  is  significantly  less  scatter  in  the  failure 
predictions  than  those  seen  in  the  boron/5505  results.  This  difference  can  be  partly 
attributed  to  the  improved  statistical  representation  of  the  experimental  failure  loads 
achieved  by  averaging  five  data  points  per  case.  MCT  and  Tsai-Wu  produced  the  same 
mean  values  for  the  E-glass/8084  uniaxial  load  cases  with  Tsai-Wu  again  producing 
significantly  more  scatter  in  failure  prediction  accuracy. 


Table  7.2  Summary  of failure  loads  for  E-glass/8084  laminates. 


Laminate 

Exp 

MCT/Exp 

Tsai-Wu/Exp 

[0/90]s  C 

-381 

1.01 

1.01 

[0/90]s  T 

374 

1.03 

1.11 

[0/90/±45]s  C 

-301 

0.95 

0.79 

[0/90/±45]s  T 

245 

1.13 

0.88 

[±45]s  C 

-110 

0.95 

0.98 

[±45]ST 

100 

0.96 

1.05 

[5]NT 

537 

1.09 

1.01 

[1  0]n  T 

281 

1.05 

1.12 

[1  5]n  T 

187 

1.04 

1.12 

[20]n  T 

140 

1.06 

1.11 

[30]n  T 

93 

1.06 

1.10 

[45]n  T 

67 

1.04 

1.04 

[60]n  T 

53 

1.06 

1.02 

Avg  = 

1.03 

1.03 

COV  = 

5.1% 

9.6% 
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7.2.3.  AS4/3501  Composite 

Failure  data  for  laminate  specimens  fabricated  from  AS4/3501  (carbon/epoxy) 
composites  were  experimentally  determined  in  previous  research  projects  conducted  by 
the  Composite  Materials  Research  Group.  Fourteen  laminates  of  different  architectures 
were  tested  under  either  uniaxial  tension  or  compression.  Stress-strain  curves  were  not 
available  for  most  of  the  AS4/3501  laminates  tested.  The  two  curves  that  were  available 
were  digitized  and  presented  here.  MCT_77  results  were  based  on  constitutive  material 
properties  backed  out  via  micromechanics  from  experimental  data  based  on  at  least  five 
samples  per  test.  Material  properties  used  in  the  analysis  are  discussed  in  Appendix  C. 
AS4/3501  lamina  have  essentially  the  same  modulus  values  in  tension  and  compression. 


Fig.  7.9  AS4/3501  [±60]s  laminate  under  uniaxial  compression. 
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The  [±60]s  laminate,  shown  in  Fig.  7.9,  initially  had  a  near  linear  load  response  to  the 
applied  compressive  load.  At  approximately  -175  MPa,  the  analytical  simulation  went 
nonlinear  due  to  significant  levels  of  shear  stress  which  appear  to  be  supported  by  the 
experimental  data.  The  MCT  curve  assumed  a  linear  stress-strain  response  above  -0.04 
shear  strain  because  shear  strain  levels  exceeded  the  stored  shear  stress  strain  curve  data. 
MCT  used  the  last  positive  tangent  shear  modulus  calculated  from  the  stored  shear  stress- 
strain  curve  to  complete  the  applied  load.  MCT  over-estimated  the  strain  response  by 
over  70%  but  estimated  the  failure  load  within  3%  of  the  experimental  data.  Tsai-Wu 
was  within  14%.  No  composite  damage  was  detected  until  compressive  failure  of  the 
fibers  precipitated  ultimate  laminate  failure. 


Fig.  7  .10  AS4/3501  [0/90] s  laminate  under  uniaxial  compression. 
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The  experimental  [0/90]s  cross-ply  laminate  data,  shown  in  Fig.  7.10,  was  slightly 
nonlinear  because  high  matrix  compressive  strength  allowed  nonlinear  matrix  behavior  to 
manifest  itself  in  the  90  plies  (see  lamina  transverse  compression  stress-strain  curve  in 
Fig.  C.4).  MCTJ77  assumed  linear  transverse  compression  behavior  but  still  predicted  a 
failure  load  within  1%  of  the  experimentally  determined  value.  The  Tsai-Wu  failure 
prediction  was  off  by  26%. 


Table  7.3  Summary  of  failure  loads  for  AS4/3 501  laminates. 


Laminate 

Exp 

MCT/Exp 

Tsai-Wu/Exp 

[±15]s  T 

786 

1.40 

2.80 

[±22]s*  T 

786 

1.01 

2.02 

[±30]s  T 

455 

1.00 

2.05 

[±45]s  T 

155 

0.95 

1.08 

[±60]s  C 

-288 

0.97 

1.14 

[±60]s  T 

74 

1.16 

0.73 

[±75]s  T 

43 

1.02 

0.77 

[90/0]3S  C 

-1074 

1.01 

1.26 

[(±45/02)2/90]s  C 

-1074 

1.03 

0.67 

[±45/0/90]2S  C 

-855 

0.91 

0.69 

[±452/03/±45]s  C 

-827 

0.94 

0.44 

[±452/03/902]s  C 

-746 

1.22 

0.97 

[±452/90/0]s  C 

-753 

0.83 

0.53 

[±452/903/02]s  C 

-599 

1.20 

1.02 

■ 

Avg  = 

1.05 

1.15 

COV  = 

14.3% 

58.8% 

A  summary  of  MCT  and  Tsai-Wu  failure  predictions  for  the  AS4/3501  laminates 
analyzed  is  presented  in  Table  7.3.  The  table  shows  MCT  produced  a  78%  improvement 
in  failure  prediction  accuracy  over  Tsai-Wu  in  average  absolute  differences.  This  was  a 
much  greater  accuracy  improvement  than  that  seen  for  boron/5505  (32%)  or  E-glass/8084 
(38%).  The  source  of  this  difference  is  not  clear,  but  it  should  be  noted  that  both  boron 
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and  E-glass  are  isotropic  fibers  whereas  AS4  carbon  is  anisotropic.  Further,  the 
boron/5505  and  E-glass  laminates  had  low  fiber  volume  fractions  (-50%),  whereas  the 
AS4/3501  laminates  had  a  significantly  higher  fiber  volume  fraction  (-67%). 

7.3.  Biaxial  Loading 

Arguments  can  be  made  that  most  composite  structures  typically  operate  under  multi- 
axial  load  states,  hence  testing  composite  laminates  under  such  conditions  would  be  the 
prudent  course  of  action.  Unfortunately,  triaxial  or  biaxial  load  tests  of  composite 
laminates  are  considerably  more  difficult  to  accomplish  than  uniaxial  load  tests.  There  is 
an  order-of-magnitude  increase  in  the  complexity  of  both  the  test  frames  which  apply  the 
load  and  the  design  and  fabrication  of  laminate  test  specimens.  As  a  result,  there  is  a 
paucity  of  multi-axial  experimental  failure  data  to  verify  analysis.  Welsh  [42]  designed 
and  fabricated  a  next-generation  triaxial  test  frame  and  experimentally  developed  biaxial 
(tension-tension  and  tension-compression)  failure  surfaces  for  AS4/3501  [0/90]s 
laminates.  The  experimentally  determined  biaxial  failure  data  along  with  two 
dimensional  failure  surfaces  developed  using  MCT  and  Tsai-Wu  failure  criteria  are 
shown  in  Fig.  7.11.  The  failure  surface  generated  by  the  MCT  failure  criterion  and  the 
experimental  data  were  in  excellent  agreement.  The  failure  surface  generated  by  the 
Tsai-Wu  failure  criteria  was  not. 
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Lack  of  accuracy  in  Tsai-Wu  predictions  can  be  attributed  to  two  fundamental  problems. 
The  first  can  be  illustrated  by  examining  two  load  paths,  shown  in  Fig.  7.12  occurring  in 
quadrant  IV  of  two-dimensional  stress  space.  On  leg  one  of  load  path  #1,  a  200  MPa 
tensile  load  is  applied  to  the  laminate  in  the  0°  direction.  Holding  the  200  MPa  load 
constant,  the  second  leg  of  the  load  path  applies  a  600  MPa  compressive  load  in  the  90° 
direction  until  compressive  fiber  failure  occurs.  Load  path  #2  is  identical  to  #1  except  the 
magnitude  of  the  first  leg  is  increase  from  200  MPa  to  300  MPa.  Tsai-Wu  now  detects 
tensile  failure  of  the  90°  plies  and  reduces  the  composite  material  stiffness  and 
subsequent  load  carry  capability  to  near  zero.  When  leg  two  of  load  path  #2  is  applied  to 
the  90°  plies,  for  the  Tsai-Wu  case  there  is  no  structural  resistance,  displacements  grow 
without  bound,  and  failure  is  declared. 

The  MCT  failure  criterion  also  detects  tensile  failure  in  the  90°  plies  during  the  first  leg 
of  path  #2,  but  determines  only  matrix  failure  has  occurred.  Composite  material 
properties  are  reduced  to  account  for  macro-damage  due  to  matrix  failure.  However,  the 
reduced  composite  properties  retain  the  undamaged  fiber’s  ability  to  sustain  load.  When 
leg  two  load  is  applied,  the  laminate  continues  to  load,  in  a  manner  almost  identical  to 
path  #1,  until  compressive  fiber  failure  occurs.  An  inability  of  the  Tsai-Wu  criterion  to 
differentiate  between  minor  (matrix)  and  major  (fiber)  constituent  damage  overly 
penalizes  any  simulation  of  a  laminate’s  load  carrying  capability. 

A  second  shortcoming  in  the  Tsai-Wu  criterion  is  found  in  the  failure  parameters  F1-F3 
and  F11-F33  (see  Section  4.1).  Their  development  results  in  an  questionable  dependence 
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of  failure  predictions  in  the  purely  compressive  quadrant  m  region  of  stress  space  on 
tensile  strength  properties  and  vise-versa  for  the  purely  tensile  quadrant  II.  Lamina 
strengths  are  clearly  dependent  on  the  tensile  or  compressive  nature  of  the  applied  load 
and  its  orientation  relative  to  the  composite’s  constituents.  The  effort  to  develop  a  single 
smooth  and  continuous  failure  surface,  while  convenient,  is  simply  too  constraining.  A 
single,  smooth  curve  cannot  adequately  account  for  the  variety  of  failure  modes  arising 
from  radical  strength  differences  between  constituents  and  significant  tension  and 
compression  strength  differences  within  each  constituent. 

Significant  differences  in  tensile  and  compressive  longitudinal  and  transverse  strengths 
(831  and  -1079  longitudinal  and  28  MPa  and  -282  MPa  transverse  respectively)  for 
AS4/3501  lamina  and  the  assumption  of  developing  a  smooth  surface  result  in  the 
skewed  Tsai-Wu  failure  surface  of  quadrants  I  and  HI.  This  effect  is  somewhat  mitigated 
for  [0/90] s  E-glass/8084  laminates,  shown  in  Fig.  7.13,  where  tensile  and  compressive 
longitudinal  and  transverse  strengths  (818  and  -759  longitudinal  and  45  MPa  and  -144 
MPa  transverse  respectively)  are  closer  together.  Interestingly,  Tsai-Wu  and  MCT 
failure  surfaces  for  a  [0/90]  s  laminate  under  combined  in-plane  shear  and  normal  loading, 
shown  in  Fig.  7.14,  are  virtually  identical. 


Failure  surfaces  for  quasi-isotropic,  [0/90/±45]s,  E-glass/8084  laminates  under  normal- 
normal  and  normal-shear  loading  are  shown  in  Fig.  7.15  and  Fig.  7.16  respectively. 
Currently  there  are  no  experimental  data  to  verify  these  surfaces. 
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90  LAMINATE  STRESS 
(MPa) 


Fig.  7.13  E-glass/8084  [0/90] s  laminate  biaxial  failure  envelopes  for  combined  0-90 
loading. 


Fig.  7.14  E-glass/8084  [0/90]  s  laminate  biaxial  failure  envelopes  for  combined  0-Shear 
loading. 
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Fig.  7.15  E-glass/8084  [0/90/±45]s  laminate  biaxial  failure  envelopes  for  combined  0-90 
loading. 


Fig.  7.16  E-glass/8084  [0/90/±45]s  laminate  biaxial  failure  envelopes  for  combined  0- 
Shear  loading. 
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SECTION  8. 

SUMMARY  AND  RECOMMENDATIONS 

8.1.  Summary 

The  purpose  of  the  research  documented  herein  was  to  develop  a  finite  element  based 
analysis  tool  for  predicting  failure  in  composite  structural  laminates.  Multicontinuum 
Theory  provided  the  theoretical  framework  with  which  to  accomplish  this  goal.  Using 
constituent  stress  information  provided  by  MCT,  a  constituent-based  stress-interactive 
failure  criterion  was  proposed.  The  failure  criterion,  functionally  identical  to  the  Tsai- 
Wu  and  Hashin  criteria,  results  in  a  continuous  but  not  smooth  failure  surface  for  a 
composite  laminate  by  recognizing  tensile  and  compressive  failure  modes  unique  to  each 
constituent. 

When  benchmarking  the  MCT  failure  criterion  against  the  Tsai-Wu  criterion  and 
experimental  data,  MCT  resulted  in  a  significant  improvement  in  average  absolute 
differences  between  predicted  and  experimentally  determined  laminate  failure  strengths 
for  all  the  uniaxial  load  cases  considered.  There  were  large  standard  deviations  in  the 
results  which  are  characteristic  of  composite  failure  phenomena.  Although  using  MCT 
resulted  in  significant  improvements  in  prediction  precision  over  Tsai-Wu  in  uniaxial 
load  cases,  their  differences  became  even  more  apparent  in  the  analysis  of  biaxial  load 
cases.  For  the  single  set  of  experimental  biaxial  data  available,  the  MCT-generated 
failure  surface  was  in  excellent  agreement  while  the  Tsai-Wu  generated  surface  was  not. 
Further,  the  Tsai-Wu  failure  surface  was  neither  consistently  conservative  nor  un- 
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conservative  when  compared  to  the  experimental  data,  which  is  unacceptable  in  structural 
design.  Since  MCT  failure  methodology  has  been  shown  to  be  successful  using  three 
dimensional  constituent  stress  fields  derived  from  uniaxial  and  biaxial  tests,  there  is  good 
reason  to  believe  that  it  will  also  be  successful  in  failure  predictions  involving  three- 
dimensional  load  cases.  Increased  failure  prediction  confidence,  provided  by  MCT, 
requires  a  few  more  material  tests  ([±0]s  to  optimize  constituent  shear  strengths)  than 
those  necessary  to  determine  the  failure  tensors  for  an  analogous  Tsai-Wu  analysis. 

* 

8.2.  Recommendations  for  Future  Research 

Experimentally  generated  stress-strain  curves  clearly  show  that  composite  behavior  is 
viscoelastic.  The  nonlinear  elastic  constituent  model  used  in  this  research  should  be 
replaced  with  a  viscoelastic  one.  Because  the  source  of  nonlinearity  in  organic 
composites  is  rooted  in  the  matrix,  a  viscoelastic  model  should  address  all  potential 
sources  of  nonlinearity  (e.g.  shear  and  transverse  compression).  Should  research 
emphasis  turn  to  metal  matrix  composites,  an  elastic-plastic  constitutive  model  should 
also  be  available. 

Multicontinuum  theory  is  not  restricted  by  the  geometry  of  underlying  microstructure.  In 
this  research,  only  lamina  based  on  hexagonal  fiber  packing  were  considered.  A  library 
of  micromechanics  models  should  be  developed  to  broaden  MCT_77’s  analytical  scope. 
The  library  could  include  microstructures  for  the  following: 

a.  Random  fiber  packing  of  unidirectional  lamina, 

b.  Woven  fabrics, 

c.  Stitched  fabrics, 
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d.  Randomly  distributed  short  (chopped)  fiber  reinforcement,  and 

e.  Randomly  distributed  particulate  reinforcement. 

Further,  Multicontinuum  Theory  could  be  extended  to  a  third  constituent  to  account  for 
interfacial  phases  between  the  matrix  and  reinforcement  or  the  warp,  weft,  and  matrix 
components  of  a  woven  fabric  composite. 


Upgrading  MCT  77  software  to  bring  it  more  inline  with  commercial  finite  element 
codes  would  be  a  worthwhile  endeavor.  Without  comment,  the  follow  upgrades  could 
enhance  MCT_77’s  analysis  capabilities 


a.  Ability  to  degenerate  elements  to  a  lower  number  of  nodes; 

b.  Ability  to  mix  element  types  within  a  single  analysis; 

c.  Parse  a  single  user  input  file  to  determine  quantities  to  be  read  and  conduct  error 
checking; 

d.  A  computationally  more  efficient  method  for  solving  matrix  equations,  such  as 
Inverse-Broyden,  should  be  incorporated  for  incremental  analyses; 

e.  General  improvements  in  the  ANSYS  to  MCT_77  translator; 

f.  Expand  element  library  (twenty-node  BRICK  and  improved  SHELL  elements); 
and 

g.  Port  MCT  77  code  to  FORTRAN  90  to  take  advantage  of  dynamic  memory 
allocation. 

Development  of  the  MCT  failure  criteria  is  not  complete.  The  criteria  could  be  expanded 
to  included  fiber  buckling  which  is  especially  significant  in  carbon/epoxy  laminates. 
Fiber-matrix  debonding  could  also  be  incorporated  as  a  failure  mode.  The  nonlinear 
regression  optimization  procedure  for  determining  constituent  shear  strengths  typically 
has  poor  convergence.  Improved  optimization  techniques  should  be  available. 


M.J.  Hinton,  P.Soden,  and  S.  Kaddour  [43]  are  conducting  a  “World  Wide  Exercise”  to 
determine  the  status  of  methods  for  predicting  the  strength  of  organic  composite 
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laminates.  Proponents  of  approximately  a  dozen  different  failure  theories  were  given 
identical  material  data  and  asked  to  predict  failure  for  a  given  set  of  laminates  under 
biaxial  load.  The  exercise  is  ongoing  and  is  scheduled  for  completion  in  1999.  The 
experimental  data  and  analytical  results  from  this  study  should  provide  invaluable  insight 
into  the  phenomena  of  composite  failure. 

The  ultimate  verification  of  MCT  or  any  other  structural  failure  criteria  is  found  in  the 
analysis  of  actual  in-service  structures.  Verification  of  this  sort  will  require  acceptance  of 
the  criterion  by  a  number  of  working  engineers,  a  time  frame  measured  in  decades,  and 
funding  in  the  tens  of  millions  of  dollars. 
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APPENDIX  A 

BORON/5505  LAMINA  PROPERTIES 

A.l.  Material  Test  Suite 

The  composite  material  properties  presented  in  this  appendix  are  from  lamina  composed 
of  boron  reinforcing  fibers  in  a  Narmco  5505  epoxy  matrix.  The  fiber  volume  fraction 
was  given  as  50%.  Experimental  stress-strain  curves  for  five  of  the  tests  making  up  the 
material  test  suite  (longitudinal  and  transverse  tension  and  compression,  and  longitudinal 
shear)  were  digitized  from  figures  published  by  Petit  and  Waddoups  [1]  and  are  presented 
in  Fig.  A.l  -  Fig.  A. 5.  Only  five  tests  of  the  material  test  suite  were  conducted  on  this 
composite.  No  out-of-plane  shear  response  was  recorded.  Only  in-plane  load  cases  were 
considered.  The  authors  note  that  the  stress-strain  curves  presented  are  based  on  a  single 
test,  hence  there  is  no  statistical  basis  for  the  results.  No  curve  fits  were  used  to  plot  the 
experimental  data,  i.e.,  the  lines  are  linear  between  the  digitized  points. 

MCT_77  results,  based  on  constitutive  material  properties  backed  out  via 
micromechanics  from  the  experimental  data,  are  also  plotted.  Boron/5505  lamina  had 
different  tension  and  compression  modulus  values.  MCT_77  currently  supports  only  one 
modulus  value  per  direction,  so  tensile  modulus  values  were  used  in  all  analyses. 

Table  A.l  presents  experimentally  determined  initial  composite  elastic  constants.  Table 
A.2  presents  in  situ  constituent  elastic  constants  backed  out  from  composite  values. 
Table  A.3  presents  curve  fit  parameters  for  the  nonlinear  composite  shear  stress-strain 
curve. 
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Table  A.l  Boron/5505  composite  elastic  constants  derived  from  experimental  data. 


Eu  (Gpa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa) 

Vi2 

FVF 

207.0 

21.0 

7.71 

7.23 

0.294 

0.50 

Table  A.2  In  situ  elastic  constants  for  boron/5505  constituents. 


Constituent 

E(GPa) 

V 

Boron 

407 

0.25 

5505  epoxy 

7.24 

0.35 

Table  A.3  Nonlinear  shear- 12  curve  fit  parameters  for  boron/5505. 


Co  (MPa) 

C,  (MPa) 

C2  (MPa) 

ai 

a2 

138.6 

-136.9 

-1.726 

-45.23 

-881.8 

Fig.  A.1  Boron/5505  lamina  under  longitudinal  tension. 
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Fig.  A.4  Boron/5505  lamina  under  transverse  compression 


Fig.  A.5  Boron/5505  lamina  under  longitudinal  shear. 
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Note  that  the  shear  stress/strain  results  were  taken  from  the  results  of  a  [±45]s  tension 
test.  The  curve  is  a  combination  of  digitized  results  taken  from  the  graphs  of  Petit  and 
Waddoups  and  tabular  stress/strain  results  published  by  Rosen  [2]  based  on  additional 
work  published  by  Petit  [3]. 


A.2.  Micromechanics  Results:  Referenced  to  Changing  Gn 

Changes  in  composite  and  matrix  elastic  constants  with  respect  to  change  in  the 
composite  shear  modulus  were  determined  by  invoking  the  constitutive  model  described 
in  section  5.3.  Although  it  is  the  composite  shear  strain  that  initiates  recalculation  of  all 
material  constants  within  the  MCT  77  code,  the  composite-constituent  relationships  were 
developed  by  varying  the  matrix  shear  modulus  from  0%  to  100%  of  its  initial  in  situ 
value  and  using  micromechanics  to  derive  new  composite  properties.  Once  a  set  of 
constituent-composite  properties  are  developed,  the  relationships  are  redefined  in  terms 
of  change  in  shear  modulus.  A  curve  for  each  elastic  constant,  normalized  with  respect  to 
the  tangent  composite  shear  modulus,  was  fitted  through  the  values  determined  at  each 
data  point.  Quadratic  equations  for  four  curves,  E22,  V12,  G23,  and  Gm,  referenced  to 
changes  in  composite  G12  are  presented  in  Fig.  A.6  through  Fig.  A.9.  Numerical  results 
are  available  in  Table  A.4. 
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y  =  -1.201E+00X2  +  3.853E+00X  +  4.270E-02 
R2  =  9.992E-01 


L6  Boron/5505:  E22  as  a  function  of  change  in  G12. 


y  =  2.374E-03X2  -  1.290E-02x  +  4.861  E-02 
R2  =  9.999E-01 


♦  Nu12 

- Quadratic  fit 


Fig.  A.7  Boron/5505:  vI2  as  a  function  of  change  in  Gn- 
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.8  Boron/5505:  G23  as  a  function  of  change  in  Gn- 


0.4000 

0.3500 

0.3000 

0.2500 

0.2000 

0.1500 

0.1000 

0.0500 

0.0000 


y  =  1.474E-02X2  +  3.327E-01X  +  1.234E-05 
R2  =  1.000E+00 


0%  20%  40%  60%  80%  100% 

G127G12 


Fig.  A.9  Boron/5505:  G  as  a  function  of  change  in  G12. 
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Table  A.4  Micromechanics  results  for  boron/5505  lamina  properties  referenced  to 
changing  Gn  shear  modulus. 


Matrix  Value  Value/G12 

Corrmosite  Value  Value/G12 

i  100%  | 

E  (GPa)  7.239 

G  (GPa)  2.681  0.3475 

Nu  0.350 

K  (GPa)  8.046 

Ell  (GPa)  206.9  26.82 

E22  (GPa)  21.00  2.722 

G12  (GPa)  7.715  1.000 

G23  (GPa)  7.229  0.9370 

Nu12  0.294  0.0381 

1  88%  | 

E  (GPa)  6.448  . 

G  (GPa)  2.359  0.3058 

Nu  0.366 

K  (GPa)  8.046 

Ell  (GPa)  206.5  26.77 

E22  (GPa)  19.24  2.494 

G12  (GPa)  6.825  0.885 

G23  (GPa)  6.489  0.841 1 

Nu12  0.301  0.0390 

I  66%  | 

E  (GPa)  4.946 

G  (GPa)  1.769  0.2294 

Nu  0.398 

K  (GPa)  8.046 

Ell  (GPa)  205.7  26.66 

E22  (GPa)  15.82  2.051 

G12  (GPa)  5.164  0.669 

G23  (GPa)  5.085  0.6591 

Nu12  0.317  0.0411 

I  44%  | 

E  (GPa)  3.374  i 

G  (GPa)  1.180  0.1529  j 

Nu  0.430 

K  (GPa)  8.046 

Ell  (GPa)  204.9  26.56 

E22  (GPa)  11.820  1.532 

G12  (GPa)  3.476  0.4506 

G23  (GPa)  3.586  0.4648 

Nu12  0.334  0.0433 

1  22%  1 

E  (GPa)  1.727 

G  (GPa)  0.590  0.0765 

Nu  0.464 

K  (GPa)  8.046 

Ell  (GPa)  204.100  26.45 

E22  (GPa)  6.901  0.8945 

G12  (GPa)  1.755  0.227 

G23  (GPa)  1.936  0.2509 

Nu12  0.353  0.0458 

1  12%  I 

E  (GPa)  0.952 

G  (GPa)  0.322  0.0417 

Nu  0.480 

K  (GPa)  8.046 

Ell  (GPa)  203.7  26.40 

E22  (GPa)  4.122  0.5343 

G12  (GPa)  0.9610  0.125 

G23  (GPa)  1.1050  0.1432 

Nu12  0.363  0.0471 

1  5%  i 

E  (GPa)  0.3999 

G  (GPa)  0.1341  0.0174 

Nu  0.492 

K  (GPa)  8.046 

Ell  (GPa)  203.4  26.36 

E22  (GPa)  1 .858  0.2408 

G12  (GPa)  0.4016  0.052 

G23  (GPa)  0.4787  0.0620 

Nu12  0.370  0.0480 

i  0.1%  i 

E  (GPa)  0.008042 

G  (GPa)  0.002681  0.0003 

Nu  0.4998 

K  (GPa)  8.046 

Ell  (GPa)  203.2  26.34 

E22  (GPa)  0.039  0.0051 

G12  (GPa)  0.0081  0.001 

G23  (GPa)  0.0099  0.0013 

Nu12  0.375  0.0486 
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A. 3.  Micromechanics  Results:  Referenced  to  Changing  G23 

Petit  and  Waddoups  did  not  publish  transverse  shear  data.  All  load  cases  published 
within  the  paper  were  in-plane  loading  only,  precluding  the  necessity  of  out-of-plane 
data.  MCT_77  analysis  conducted  on  these  load  cases  assumed  in-plane  and  out-of-plane 
shear  stress-strain  curves  to  be  identical  in  shape.  Ultimate  constituent  shear  strengths, 
both  in-  and  out-of-plane  are  determined  analytically  using  the  optimization  procedure 
outlined  in  section  5.4. 


A. 4.  Failure  Data:  Elastic  Constants 

Failure  of  a  constituent  is  assumed  absolute  without  directionality.  Moduli  for  the  failed 
constituent  are  set  at  near  zero  values,  rather  than  zero,  to  avoid  numerical  difficulties 
within  the  MCT_77  code.  As  shown  in  Table  A.5,  matrix  moduli  are  reduced  to  1%  of 
original  value.  Failed  fiber  moduli,  which  are  initially  one  to  two  orders  of  magnitude 
larger  than  those  of  the  matrix,  are  reduced  such  that  their  magnitudes  are  of  the  same 
order  as  that  of  failed  matrix  to  avoid  any  artificial  stiffness  in  failed  regions.  Damaged 
composite  properties  that  reflect  failed  constituents,  shown  in  Table  A.6  and  Table  A.7, 
are  calculated  via  the  micromechanics  model.  A  complete  input  file  of  material 
properties  for  a  MCT  77  failure  analysis  consist  of  seven  materials:  undamaged 
composite,  intact  fiber,  intact  matrix,  fiber-failed  composite,  matrix-failed  composite, 
failed  fiber,  and  failed  matrix.  A  fiber-failed  composite  assumes  simultaneous  matrix 
failure  also. 
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Table  A.5  Elastic  constants  for  boron/5505  constituents  assuming  failure. 


Constituent 

E  (GPa) 

%  of  Original  Value 

V 

Boron 

0.0407 

0.01 

0.25 

5505  epoxy 

0.0724 

1.0 

0.35 

Table  A.6  Boron/5505  composite  elastic  constants  with  failed  fiber/matrix. 


En  (GPa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa) 

V12 

FVF 

0.0567 

0.0540 

0.0210 

0.0208 

0.308 

0.50 

Table  A.7  Boron/5505  composite  elastic  constants  with  failed  matrix. 


E„  (GPa) 

E22  (GPa) 

G12  (GPa) 

G23  (GPa) 

V12 

FVF 

203.0 

0.220 

0.0805 

0.00752 

0.293 

0.50 

A.5.  Failure  Data:  Ultimate  Strengths 

The  ultimate  strengths  of  the  composite  are  taken  from  material  test  suit  results.  Ultimate 
strengths  of  fiber  and  matrix  are  backed  out  of  composite  results  via  a  MCT_77  analysis 
using  material  constants  listed  in  Table  A.l  through  Table  A.3  and  the  composite 
strengths  listed  in  Table  A.8.  Determination  of  constituent  shear  strengths  are  further 
refined  by  nonlinear  regression  optimization  of  experimental  results  from  [+0]  tests  as 
discussed  in  section  5.4.  Composite  ultimate  strengths  are  not  used  directly  in  the  failure 
criterion  implemented  in  MCT_77. 
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Table  A.8  Boron/5505  composite  ultimate  strengths. 


+S„  (MPa) 

~S„  (MPa) 

~S22  (MPa) 

*^i2  (MPa) 

S23  (MPa) 

1379.0 

2864.0 

87.6 

262.0 

125.5 

95.52* 

*S23  determined  analytically  by  nonlinear  regression. 


Table  A.9  In  situ  boron  fiber  strengths. 


+Slf  (MPa) 

S"'  (MPa) 

Si V  (MPa) 

2709 

-6454 

180 

Table  A.10  In  situ  5505  matrix  strengths. 


+  o  22m 
°22 

(MPa) 

-  o  22  m 
^22 

(MPa) 

+  ri22m 
°33 

(MPa) 

-  r*  22m 

°33 

(MPa) 

r»12m 

°12 

(MPa) 

r»23m 

^23 

(MPa) 

70.9 

-212 

7.53 

-22.5 

83.2 

64.2 

Constituent  ultimate  strengths  are  listed  in  Table  A.9  and  Table  A.  10.  The  subscript  is 
the  component  of  ultimate  strength  and  the  superscript  is  the  direction  of  ultimate  applied 
load,  e.g.,  +Sf  is  the  strength  (max  stress)  in  the  33  direction  when  an  ultimate  22  tensile 

load  is  applied.  All  constituent  ultimate  strengths  listed  are  required  input  to  the  failure 
criterion  implemented  in  MCT_77. 


A.  6.  Uniaxial  Laminate  Tests 

Experimental  stress-strain  curves  for  laminate  specimens  fabricated  from  boron/5505 
(boron/epoxy),  shown  in  Fig.  A.  10  -  Fig.  A.24,  were  digitized  from  figures  published  by 
Petit  and  Waddoups.  The  authors  note  that  the  stress-strain  curves  presented  are  based  on 
a  single  test,  hence  there  is  no  statistical  basis  for  the  results.  No  curve  fits  were  used  to 
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plot  the  experimental  data  against  MCT_77  results,  i.e.,  the  lines  are  linear  between 
digitized  points. 


Analytical  results  based  on  a  non-finite  element  program  developed  by  Petit  and 
Waddoups  are  also  presented.  Their  program  calculated  laminate  stress  at  a  point  based 
on  tangent  moduli  from  five  lamina  curves  (longitudinal  and  transverse  lamina  tension 
and  compression,  and  in-plane  shear)  at  each  step  of  an  incrementally  applied  load.  Their 
approach  is  remarkably  accurate,  especially  considering  the  work  was  conducted  in  the 
late  1960’s. 


Fig.  A.10  Boron/5505  [0/90] s  laminate  under  uniaxial  tension. 
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ig.  A.23  Boron/5505  [±60]s  laminate  under  uniaxial  compression. 


Fig.  A.24  Boron/5505  [(0)i/±45]s  laminate  under  uniaxial  compression. 
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APPENDIX  B 

E-GLASS/8084  LAMINA  PROPERTIES 

B.l.  Material  Test  Suite 

The  composite  material  properties  presented  in  this  appendix  are  from  lamina  composed 
of  E-glass  reinforcing  fibers  in  a  Dow  Derakane  8084  rubber-toughened  vinylester  matrix 
with  a  proprietary  modification  by  Seemann  Composites.  Average  lamina  fiber  volume 
fraction  is  51%.  Stress-strain  curves  for  six  tests  making  up  the  material  test  suite 
(longitudinal  and  transverse  tension  and  compression,  and  longitudinal  and  transverse 
shear)  are  presented  in  Fig.  B.l  through  Fig.  B.6.  These  stress-strain  curves  were 
experimentally  determined  by  the  University  of  Wyoming’s  Composite  Materials 
Research  Group  from  unidirectional  laminates  provided  by  Seemann  Composites,  Inc. 
All  laminates  were  fabricated  using  a  proprietary  vacuum  assisted  resin  transfer  mold 
process  known  as  Seemann  Composites  Resin  Infusion  Process  (SCRIMP  ). 

The  E-glass/8084  results  presented  are  based  on  at  least  five  samples  per  test.  Table  B.l 
presents  experimentally  determined  initial  composite  elastic  constants.  Table  B.2 
presents  in  situ  constituent  elastic  constants  backed  out  from  composite  values.  Table 
B.3  and  Table  B.4  presents  curve  fit  parameters  for  composite  nonlinear  shear  stress- 
•  strain  curves. 


Table  B.l  E-glass/8084  composite  elastic  constants  derived  from  experimental  data. 


Eh  (GPa) 

E22  (GPa) 

Gi2  (GPa) 

G23  (GPa) 

v)2 

FVF 

E-glass/8084 

38.5 

11.7 

4.74 

4.34 

0.273 

0.51 
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Table  B.2  In  situ  elastic  constants  for  E-glass/8084  constituents. 


Constituent 

E  (GPa) 

V 

E-glass 

71.0 

0.260 

8084  vinylester 

4.66 

0.292 

Table  B.3  Nonlinear  shear-12  curve  fit  parameters  for  E-glass/8084. 


Co  (MPa) 

C,  (MPa) 

C2  (MPa) 

ai 

a2 

60.68 

-17.68 

-43.00 

-164.7 

-42.50 

Table  B.4  Nonlinear  shear-23  curve  fit  parameters  for  E-glass/8084. 


Co  (MPa) 

Ci  (MPa) 

C2  (MPa) 

ai 

a2 

•  50.66 

-1320.0 

1269.0 

-77.37 

-77.04 

Fig.  B.l  E-glass/8084  lamina  under  longitudinal  tension. 
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Fig.  B.6  E-glass/8084  lamina  under  transverse  shear. 

B.2.  Micromechanics  Results:  Referenced  to  Changing  Gu 

Changes  in  composite  and  matrix  elastic  constants  with  respect  to  change  in  the 
composite  shear  modulus  were  determined  by  invoking  the  constitutive  model  described 
in  section  5.3.  Although  it  is  the  composite  shear  strain  that  initiates  recalculation  of  all 
material  constants  within  the  MCT  77  code,  the  composite-constituent  relationships  were 
developed  by  varying  the  matrix  shear  modulus  from  0%  to  100%  of  its  initial  in  situ 
value  and  using  micromechanics  to  derive  new  composite  properties.  Once  a  set  of 
constituent-composite  properties  are  developed,  the  relationships  are  redefined  in  terms 
of  change  in  shear  modulus.  A  curve,  for  each  elastic  constant  normalized  with  respect  to 
the  tangent  shear  modulus,  was  fitted  through  the  values  determined  at  each  data  point. 
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Quadratic  equations  for  four  curves,  E22,  vi2,  G23,  and  Gm,  referenced  to  changes  in 
composite  G,2  are  presented  in  Fig.  B.7  through  Fig.  B.10.  Numerical  results  are 
available  in  Table  B.5. 


Fig.  B.7  E-glass/8084:  E22  as  a  function  of  change  in  G12. 


Fig.  B.8  E-glass/8084:  v,2  as  a  function  of  change  in  G!2. 
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Fig.  B.9  E-glass/8084:  G23  as  a  function  of  change  in  G  12- 


Fig.  B.10  E-glass/8084:  G  as  a  function  of  change  in  Gj 
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B.3.  Micromechanics  Results:  Referenced  to  Changing  G23 

Quadratic  equations  for  four  curves,  E22,  V12,  G12,  and  Gm,  referenced  to  changes  in 
composite  G23  are  presented  in  Fig.  B.ll  through  Fig.  B.14.  Numerical  results  are 
available  in  Table  B.6. 


Fig.  B.ll  E-glass/8084:  E22  as  a  function  of  change  in  G23. 


Fig.  B.12  E-glass/8084:  V\2  as  a  function  of  change  in  G23. 
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Table  B.6  Micromechanics  results  for  E-glass/8084  lamina  properties  referenced  to 
changing  G23  shear  modulus. 

Matrix  Values  Value/G23  [ 
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B.4.  Failure  Data:  Elastic  Constants 

Failure  of  a  constituent  is  assumed  absolute  without  directionality.  Moduli  for  failed 
constituents  are  set  at  near  zero  values,  rather  than  zero,  to  avoid  numerical  difficulties 
within  the  MCT_77  code.  As  shown  in  Table  B.7,  matrix  moduli  are  reduced  to  1%  of 
original  value.  Failed  fiber  moduli,  which  are  initially  one  to  two  orders  of  magnitude 
larger  than  those  of  the  matrix,  are  reduced  such  that  their  magnitudes  are  of  the  same 
order  as  that  of  failed  matrix  to  avoid  any  artificial  stiffness  in  failed  regions.  Damaged 
composite  properties  that  reflect  failed  constituents,  shown  in  Table  B.8  and  Table  B.9, 
are  calculated  via  the  micromechanics  model.  A  complete  input  file  of  material 
properties  for  a  MCT_77  failure  analysis  consist  of  seven  materials:  undamaged 
composite,  intact  fiber,  intact  matrix,  fiber-failed  composite,  matrix-failed  composite, 
failed  fiber,  and  failed  matrix.  A  fiber-failed  composite  assumes  simultaneous  matrix 

failure  also. 


Table  B.7  Elastic  constants  for  E-glass/8084  constituents  assuming  failure. 


Constituent 

E  (Gpa) 

%  of  Original  Value 

V 

E-glass 

0.071 

0.1 

0.26 

8084  vinylester 

r  0.047 

1.0 

0.292 

Table  B.8  E-glass/8084  composite  elastic  constants  with  failed  fiber/matrix. 


En  (GPa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa) 

V12 

FVF 

E-glass/8084 

0.0590 

0.0575  . 

0.0226 

0.224 

0.275 

0.51 

Table  B.9  E-glass/8084  composite  elastic  constants  with  failed  matrix. 


En  (GPa) 

E22  (GPa) 

G12  (GPa) 

G23  (GPa) 

V,2 

FVF 

E-glass/8084 

36.2  1 

0.135 

0.0555 

0.0495 

0.279 

0.51 
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B.5.  Failure  Data:  Ultimate  Strengths 

Ultimate  strengths  of  the  composite  are  taken  from  material  test  suit  results.  Ultimate 
strengths  of  the  fiber  and  matrix  are  backed  out  of  composite  results  via  a  MCT_77 
analysis  using  material  constants  listed  in  Table  B.l  through  Table  B.3  and  the  composite 
strengths  listed  in  Table  B.10.  Determination  of  constituent  shear  strengths  are  further 
refined  by  nonlinear  regression  optimization  of  experimental  results  from  [±0]  tests  as 
discussed  in  section  5.4.  Composite  ultimate  strengths  are  not  used  directly  in  by  the 
failure  criterion  implemented  in  MCT_77. 


Table  B.10  E-glass/8084  composite  ultimate  strengths. 


+5U  (MPa) 

~Sn  (MPa) 

+S22  (MPa) 

~S22  (MPa) 

£,2  (MPa) 

S23  (MPa) 

817.5 

-759.7 

45.26 

-144.3.0 

48.52 

Table  B.ll  In  situ  E-glass  fiber  strengths. 


Xv  (MPa) 

~S"f  (MPa) 

S\lf  (MPa) 

1507 

-1399 

120 

Table  B.l  2  In  situ  8084  matrix  strengths. 


ol2m 

0]2 

(MPa) 

r»23m 

°23 

(MPa) 

37.10 

34.39 

25.21 

The  constituent  ultimate  strengths  are  listed  in  Table  B.ll  and  Table  B.12.  The  subscript 
is  the  component  of  ultimate  strength  and  the  superscript  is  the  direction  of  ultimate 
applied  load,  e.g.,  +Sf  is  the  strength  (max  stress)  in  the  33  direction  when  an  ultimate 
22  tensile  load  is  applied.  All  constituent  ultimate  strengths  listed  are  required  input  to 
the  failure  criterion  implemented  in  MCT  77. 
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B.6.  Uniaxial  Laminate  Tests 

Stress-strain  curves  for  laminate  specimens  fabricated  from  E-glass/8084 
(fiberglass/vinylester)  composites,  shown  in  Fig.  B. 14-Fig.  B.21,  were  experimentally 
determined  by  the  University  of  Wyoming’s  Composite  Materials  Research  Group  from 
laminates  provided  by  Seemann  Composites,  Inc.  All  laminates  were  fabricated  using  a 
proprietary  vacuum  assisted  resin  transfer  mold  process  known  as  Seemann  Composites 

Resin  Infusion  Process  (SCRIMP0). 


Fig.  B.15  Off-angle  E-glass/8084  lamina  under  uniaxial  tension. 
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Fig.  B.20  E-glass/8084  [0/90] s  laminate  under  uniaxial  compression. 


Fig.  B.21  E-glass/8084  [0/90/+45] s  laminate  under  uniaxial  compression. 
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APPENDIX  C 

CARBON/3501  LAMINA  PROPERTIES 

C.l.  Material  Test  Suite 

Composite  material  properties  presented  in  this  appendix  are  from  lamina  composed  of 
AS4  carbon  fibers  in  a  Hercules  3501-6  epoxy  matrix.  Average  lamina  fiber  volume 
fraction  ranges  from  64%  to  67%  depending  on  the  data  source.  Stress  strain  curves  for 
the  six  tests  making  up  the  material  test  suite  (longitudinal  and  transverse  tension  and 
compression,  and  longitudinal  and  transverse  shear)  are  presented  in  Fig.  C.l  through 
Fig.  C.6.  These  stress-strain  curves  were  experimentally  determined  by  the  University  of 
Wyoming’s  Composite  Materials  Research  Group  [1,2, 3, 4, 5]  and  the  Marine  Composites 
Branch,  Structures  and  Composite  Department  of  the  Naval  Surface  Warfare  Center’s 
Carderock  Division  [6,7].  Transverse  tensile  data  was  assumed  completely  linear  and 
was  taken  from  published  values  [8]. 

Most  of  the  AS4/3501  results  presented  are  based  on  least  five  samples  per  test.  Table 
C.l  presents  the  experimentally  determined  initial  elastic  constants  for  the  composite. 
Table  C.2  and  Table  C.3  presents  constituent  elastic  constants  backed  out  from  the 
composite  values.  Table  C. 4  and  Table  C.5  presents  curve  fit  parameters  for  the 
nonlinear  shear  curves. 


Table  C.l  AS4/3501  composite  elastic  constants  derived  from  experimental  data. 


Ei,  (GPa) 

E22  (GPa) 

Gn  (GPa) 

G23  (Gpa) 

v,2 

FVF 

AS4/3501 

134.0  “1 

9.14 

7.29 

3.16 

0.261 

0.66 

C-l 


NSWCCD-65-TR-1999/1 5 


Table  C.2  In  situ  elastic  constants  for  AS4  fiber . 


Constituent 

B,,  (GPa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa) 

Vl2 

V23 

AS4  carbon 

201.0 

13.5 

95.0 

4.90 

0.22 

0.25 

Table  C.3  In  situ  elastic  constants  for  3501  matrix. 


Constituent 

E  (Gpa) 

V 

3501-6  epoxy 

4.30 

0.35 

Table  C.4  Nonlinear  shear-12  curve  fit  parameters  for  AS4/3 501. 


Co  (MPa) 

Ci  (MPa) 

C2  (MPa) 

ai 

a2 

95.17 

-4189. 

1  4093. 

-32.35 

-31.32 

Table 

C.5  Nonlinear  shear-23  curve  fit  parameters  for  AS4/3501 . 

_ _ _ — - - - - - 1 

Co  (MPa) 

C,  (MPa) 

C2  (MPa) 

ai 

a2 

360.7 

-4781. 

4421, 

-0.2796 

1"  4.131 

C-2 
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Fig.  C.2  AS4/3501  lamina  under  longitudinal  compression 


Fig.  C.3  AS4/3501  lamina  under  transverse  tension. 


NSWCCD-65-TR-1 999/1 5 


Fig.  C.4  AS4/3501  lamina  under  transverse  compression. 


Fig.  C.5  AS4/3501  lamina  under  longitudinal  shear. 


NSWCCD-65-TR-1 999/1 5 


Fig.  C.6  AS4/3501  lamina  under  transverse  shear. 


C.2.  Micromechanics  Results:  Referenced  to  Changing  Gn 

Changes  in  composite  and  matrix  elastic  constants  with  respect  to  change  in  the 
composite  shear  modulus  were  determined  by  invoking  the  constitutive  model  described 
in  section  5.3.  Although  it  is  the  composite  shear  strain  that  initiates  recalculation  of  all 
material  constants  within  the  MCT_77  code,  the  composite-constituent  relationships  were 
developed  by  varying  the  matrix  shear  modulus  from  0%  to  100%  of  its  initial  in  situ 
value  and  using  micromechanics  to  derive  new  composite  properties.  Once  a  set  of 
constituent-composite  properties  are  developed,  the  relationships  are  redefined  in  terms 
of  change  in  shear  modulus.  A  curve  for  each  elastic  constant,  normalized  with  respect  to 
the  tangent  shear  modulus,  was  fitted  through  the  values  determined  at  each  data  point. 
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Quadratic  equations  for  four  curves,  E22.  vl2,  G23,  and  Gm,  referenced  to  changes  in 
composite  Gl2  are  presented  in  Fig.  C.7  through  Fig.  C.10.  Numerical  results  are 

available  in  Table  C.6. 


C-6 


NSWCCD-65-TR-1999/15 


i  - _  _ n  -  a  go nrjr?y2 

.012E-01X  +  3.642E-05 
.000E+00 

0.45UU  - 

- 

0.4000  - 

0.3500  - 

0.3000  - 

♦  G 

?!  0.2500  - 

P  A  OAAn 

—Quadratic  tit 

0  0.2000 

0.1500  - 

" 

0.1000  - 

0.0500  - 

0.0000  < 

0%  20%  40%  60%  80%  100% 

G127G12 

Fig.  C.10  AS4/3501:  Gas  a  function  of  change  in  Gn- 
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Table  C.6  Micromechanics  results  for  AS4/3501  lamina  properties  referenced  to 
changing  G12  shear  modulus. 


Matrix  Values  Value/G12 

Comoosite  Values  Value/G12 

|  100%  | 

E  (GPa)  4.295 

G(GPa)  1.591  0.2182 

Nu  0.350 

K  (GPa)  4.772 

Ell  (GPa)  134.00  18.38 

E22(GPa)  9.140  1.254 

G12  (GPa)  7.291  1.000 

G23  (GPa)  3.163  0.4338 

Nu12  0.261  0.0358 

f  80%  I 

E  (GPa)  3.506 

G(GPa)  1-273  0.1745 

Nu  0.378 

K  (GPa)  4.772 

Ell  (GPa)  133.80  18.35 

E22  (GPa)  8.497  1.165 

G12  (GPa)  5.919  0.8118 

G23  (GPa)  2.869  0.3935 

Nu12  0.271  0.0372 

|  60%  \ 

E  (GPa)  2.684 

G  (GPa)  0.954  0.1309 

Nu  0.406 

K  (GPa)  4.772 

Ell  (GPa)  133.50  18.31 

E22  (GPa)  7.685  1.054 

G12  (GPa)  4.509  0.6184 

G23  (GPa)  2.519  0.3455 

Nu12  0.280  0.0384 

1  40%  |  ""  1 

E  (GPa)  1.828 

G  (GPa)  0.636  0.0873 

Nu  0.436 

K  (GPa)  4.772 

Ell  (GPa)  133.2  18.27 

E22  (GPa)  6.607  0.906 

G12  (GPa)  3.053  0.4187 

G23  (GPa)  2.078  0.2850 

Nu12  0.291  0.0399 

5  20%  | 

E  (GPa)  0.9337 

G  (GPa)  0.318  0.0436 

Nu  0.467 

K  (GPa)  4.772 

Ell  (GPa)  132.9  18.23 

E22  (GPa)  4.895  0.6714 

G12  (GPa)  1.551  0.2127 

G23  (GPa)  1.445  0.1982 

Nu12  0.302  0.0414 

l  10%  1 

E  (GPa)  0.4720 

G(GPa)  0.1591  0.0218 

Nu  0.484 

K  (GPa)  4.772 

Ell  (GPa)  132.7  18.20 

E22  (GPa)  3.397  0.4659 

G12  (GPa)  0.7811  0.1071 

G23  (GPa)  0.9476  0.1300 

Nu12  0.309  0.0424 

1  5%  | 

E  (GPa)  0.2373 

G  (GPa)  0.0795  0.0109 

Nu  0.492 

K  (GPa)  4.772 

Ell  (GPa)  132.70  18.20 

£22  (GPa)  2.144  0.2941 

G12  (GPa)  0.3922  0.0538 

G23  (GPa)  0.5732  0.0786 

Nu12  0.312  0.0428 

1  0.1000%  | 

E  (GPa)  0.004772 

G  (GPa)  0.0016  0.0002 

Nu  0.4998 

K  (GPa)  4.772 

Ell  (GPa)  132.60  18.19 

£22  (GPa)  0.05909  0.0081 

G12  (GPa)  0.00788  0.0011 

G23  (GPa)  0.0148  0.0020 

Nu12  0.315  0.0432 
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C.3.  Micromechanics  Results:  Referenced  to  Changing  G2s 


Quadratic  equations  for  four  curves,  E22,  v!2,  G12,  and  Gm,  referenced  to  changes  in 
composite  G23  are  presented  in  Fig.  C.ll  through  Fig.  C.14.  Numerical  results  are 
available  in  Table  C.7. 
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Fig.  C.13  AS4/3501:  Gn  as  a  function  of  change  in  G23. 


Fig.  C.14  AS4/3501:  Gas  a  function  of  change  in  G23. 
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Table  C.7  Micromechanics  results  for  AS4/3501  lamina  properties  referenced  to 
changing  G23  shear  modulus. 


Matrix  Values  Value/G23 

ComDOsite  Values  Value/G23 

|  100%  “1 

E  (GPa)  4.295  1 

G  (GPa)  1.591  0.5029  1 

Nu  0.350  1 

K  (GPa)  4.772  1 

Ell  (GPa)  134.0  42.3648 

=22  (GPa)  9.140  2.8897 

312  (GPa)  7.291  2.3051 

323  (GPa)  3.163  1.0000 

Nu12  0.261  0.0825 

|  80%  1 

E  (GPa)  3.506 

G  (GPa)  1.273  0.4023 

Nu  0.378 

K  (GPa)  4.772 

Ell  (GPa)  133.8  42.3016 

E22  (GPa)  8.497  2.6864 

G12  (GPa)  5.919  1.8713 

G23  (GPa)  2.869  0.9071 

Nu12  0.271  0.0857 

!  60%  i 

E  (GPa)  2.684 

G  (GPa)  0.954  0.3018 

Nu  0.406 

K  (GPa)  4.772 

Ell  (GPa)  133.5  42.2068 

E22  (GPa)  7.685  2.4297 

G12  (GPa)  4.509  1.4255 

G23  (GPa)  2.519  0.7964 

Nu12  0.280  0.0885 

1  40%  1 

E  (GPa)  1.828 

G  (GPa)  0.636  0.2012 

Nu  0.436 

K  (GPa)  4.772 

Ell  (GPa)  133.2  42.1119 

E22  (GPa)  6.607  2.0888 

G12  (GPa)  3.053  0.9652 

G23  (GPa)  2.078  0.6570 

Nu12  0.291  0.0920 

1  20%  I 

E  (GPa)  0.934 

G  (GPa)  0.318  0.1006 

Nu  0.467 

K  (GPa)  4.772 

Ell  (GPa)  132.9  42.0171 

E22  (GPa)  4.895  1.5476 

G12  (GPa)  1.551  0.4904 

G23  (GPa)  1.445  0.4568 

Nu12  0.302  0.0955 

1  10% — J 

E  (GPa)  0.4720 

G  (GPa)  0.1591  0.0503 

Nu  0.484 

K  (GPa)  4.772 

Ell  (GPa)  132.7  41.9538 

E22  (GPa)  3.397  1.0740 

G12  (GPa)  0.781  0.2469 

G23  (GPa)  6.948  0.2996 

Nu12  0.309  0.0977 

1  5%  I 

E  (GPa)  0.2373 

G  (GPa)  0.0795  0.0251 

Nu  0.492 

K  (GPa)  4.772 

Ell  (GPa)  132.7  41.9538 

E22  (GPa)  2.144  0.6778 

G12  (GPa)  0.392  0.1240 

G23  (GPa)  0.573  0.1812 

Nu12  0.312  0.0986 

1  0.10%  I 

E  (GPa)  0.004772 

G  (GPa)  0.0016  0.0005 

Nu  0.4998 

K  (GPa)  4.772 

Ell  (GPa)  132.6  41.9222 

E22  (GPa)  0.059  0.0187 

G12  (GPa)  0.008  0.0025 

G23  (GPa)  0.015  0.0047 

Nu12  0.315  0.0996 

C-ll 
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C.4.  Failure  Data:  Elastic  Constants 

Failure  of  a  constituent  is  assume'd  absolute  without  directionality.  Moduli  for  failed 
constituents  are  set  to  near  zero  values,  rather  than  zero,  to  avoid  numerical  difficulties 
within  the  MCT_77  code.  As  shown  in  Table  C.8,  matrix  moduli  are  reduced  to  1%  of 
the  original  value.  Failed  fiber  moduli,  which  are  initially  one  to  two  orders  of 
magnitude  larger  than  those  of  the  matrix,  are  reduced  such  that  their  magnitudes  are  of 
the  same  order  as  that  of  failed  matrix  to  avoid  any  artificial  stiffness  in  failed  regions. 
Damaged  composite  properties  that  reflect  failed  constituents,  shown  in  Table  C.9  and 
Table  C.10,  are  calculated  via  the  micromechanics  model.  A  complete  input  file  of 
material  properties  for  a  MCT_77  failure  analysis  consist  of  seven  materials;  undamaged 
composite,  intact  fiber,  intact  matrix,  fiber  failed  composite,  matrix  failed  composite, 
failed  fiber,  failed  matrix.  A  fiber  failed  composite  assumes  simultaneous  matrix  failure 

also. 


Table  C.8  Elastic  constants  for  AS4/3501  constituents  assuming  failure. 


Elastic  Constant 

Value  (GPa) 

%  of  Original  Value 

K 

0.0201 

0.01% 

Ef 

E-22 

0.0135 

0.1% 

- g' - 1 

0.0950 

0.1% 

gL 

0.0490 

1.0% 

Em 

0.0430 

’  1.0% 

Table  C.9  AS4/3501  composite  elastic  constants  with  failed  fiber/matrix. 


En  (GPa) 

E22  (GPa) 

G,2  (GPa) 

G23  (GPa). 

vi2 

0.0280 

0.0201 

0.0443 

0.0309 

0.286 

Table  C.10  AS4/3501  composite  elastic  constants  with  failed  matrix. 


Eh  (GPa) 

E22  (GPa) 

Gi2  (GPa) 

G23  (GPa) 

Vl2 

132.6 

2.205 

0.0787 

0.0784 

0.257 
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C.5.  Failure  Data:  Ultimate  Strengths 

The  ultimate  strengths  of  the  composite  are  taken  from  material  test  suit  results.  Ultimate 
strengths  of  the  fiber  and  matrix  are  backed  out  of  the  composite  results  via  a  MCT_77 
analysis  using  the  material  constants  listed  in  Table  C.l  through  Table  C.3  and  the 
composite  strengths  listed  in  Table  C.ll.  Determination  of  constituent  shear  strengths 
are  further  refined  by  nonlinear  regression  optimization  of  experimental  results  from  [±0] 
tests  as  discussed  in  section  5.4.  Composite  ultimate  strengths  are  not  used  directly  in  by 
the  failure  criterion  implemented  in  MCT_77. 


Table  C.ll  AS4/3501  composite  ultimate  strengths. 


+5U  (MPa) 

'S„  (MPa) 

*S22  (MPa) 

-S22  (MPa) 

5,2  (MPa) 

523  (MPa) 

1335 

-1992 

28 

-282 

115  1 

33 

Table  C.12  In  situ  AS4  carbon  fiber  strengths. 


+Sjlf  (MPa) 

-s;y  (MPa) 

S,12/  (MPa) 

2000. 

-2984.  1 

170. 

Table  C.13  In  situ  3501  matrix  strengths. 


+  o  22  m 
°22 

(MPa) 

-  r»22m 
^22 

(MPa) 

+  r»22m 
°33 

(MPa) 

—  ft  22m 
*^33 

(MPa) 

r»12m 

°12 

(MPa) 

o  23m 
°23 

(MPa) 

22.6 

-228 

-3.63 

36.7 

45.1 

25.0 

The  constituent  ultimate  strengths  are  listed  in  Table  C.12  and  Table  C.13.  The  subscript 
is  the  component  of  ultimate  strength  and  the  superscript  is  the  direction  of  ultimate 
applied  load,  e.g.,  +S3232  is  the  strength  (max  stress)  in  the  33  direction  when  an  ultimate 
22  tensile  load  is  applied.  All  constituent  ultimate  strengths  listed  are  required  input  to 
the  failure  criterion  implemented  in  MCT_77. 
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APPENDIX  D 

USING  EXCEL  FOR  CURVE  FITTING 

Two  types  of  curve  fits  are  used  by  the  material  model  within  the  MCT_77  finite  element 

code.  The  first  determines  an  equation  defining  nonlinear  shear  stress-strain 

relationships.  The  second  defines  relationships  between  composite  and  constituent 

properties  due  a  varying  composite  shear  modulus. 

* 

D.l.  Curve  Fitting  Nonlinear  Shear  Data 

Three  terms  of  an  exponential  series  are  used  to  determine  nonlinear  composite  shear 
stress-strain  relationships.  The  equation  has  the  form 

r  =  C0+C,  „<«•>  +  C!C<-'>_  (D.i) 

where  C,  and  a,  are  curve  fit  constants,  j is  engineering  shear  strain,  and  ns  shear  stress. 
Equation  (D.l)  has  five  curve  fit  parameters  to  be  determined  for  a  best  fit  of  shear  stress- 
strain  data  points. 

Microsoft®  Excel,  in  conjunction  with  the  Solver  analysis  option,  offers  a  convenient 
environment  for  implementing  a  nonlinear  regression  approach  to  solve  for  the  five  curve 
fit  constants.  Microsoft®  Excel  Solver  is  an  Add-in  analysis  package  included  with  the 
standard  spreadsheet  routines.  If  Solver  doesn’t  appear  under  the  Tools  menu  on  the 
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toolbar,  then  it  must  be  installed  by  checking  the  Solver  option  box  under  Tools|Add- 
Ins|Available  Add-Ins. 


The  method  begins  by  placing  stress-strain  data  in  adjacent  columns  on  the  spreadsheet. 
Five  cells  are  designated  to  hold  the  five  constants  to  be  determined.  For  clarity,  labels 
identifying  the  constants  are  placed  next  to  appropriate  cells.  It  is  also  very  advantageous 
to  use  the  Insert|Name|Define  command  sequence  from  the  toolbar  to  name  each  cell  for 
the  appropriate  constant.  Subsequent  equations  using  the  constants  can  be  written  in  a 
more  natural  and  easier  to  read  manner.  See  Fig.  D.l  for  an  example  spreadsheet. 


Fig.  D.l  Nonlinear  spreadsheet  setup. 


Using  equation  (D.l),  a  formula  for  calculating  predicted  stress  values  is  typed  into  the 
first  cell  of  the  column  adjacent  to  the  stress-strain  data.  In  the  formula,  cell  reference 
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notation  ,  e.g.,  “A8’\  is  used  in  place  of  the  actual  strain  value  from  that  data  point.  For 
example,  the  formula  might  look  like 

=  c_0  +  c_l  *  EXP(a_l  *  A10)  +  c_2  *  EXP(a_2  *  A10) .  (D.2) 

The  formula  defined  for  the  first  cell  can  then  be  copied  and  pasted  to  remaining  cells. 

The  next  adjacent  column  will  contain  the  square-of-the-differences  between  predicted 
and  actual  stress  values.  Using  cell  reference  notation,  the  formula  might  look  like 


=  (B10-C10)A2 .  (D-3) 

Again,  the  formula  can  be  written  for  the  first  cell  in  the  column  and  copied  to  the 
remaining  ones. 

At  the  bottom  of  the  column  containing  the  square-of-the-differences,  the 
Insert|Function|Math&Trig|Sum  command  sequence  is  used  to  sum  the  cells  in  the 
column.  This  cell,  known  as  the  residual,  becomes  the  target  of  subsequent  Solver 
minimization  operations. 

It  is  critical  that  the  initial  shear  modulus,  defined  by  the  derivative  of  the  stress-strain 
curve,  is  equal  to  the  shear  modulus  calculated  by  the  micromechanics  model.  Designate 
one  cell  on  the  spreadsheet  to  hold  the  shear  modulus  as  calculated  by  micromechanics. 
In  another  cell,  insert  the  derivative  of  the  equation  (D.l)  evaluated  at  y=0.  Equivalence 
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of  the  shear  modulus  determined  by  micromechanics  and  the  curve  equation  first 
derivative  will  form  a  constraint  on  the  Solver  optimization. 


Results  of  a  nonlinear  analysis  regression  are  highly  dependent  on  starting  values. 
Generally,  nonlinear  composite  shear  curves  have  the  same  basic  shape,  hence,  curve  fit 
constants  that  work  for  one  data  set  provides  excellent  starting  values  for  another.  Trial 

starting  values  are  provided  in  Table  D.l. 


Table  D.l  Trial  starting  values  for  nonlinear  regression  analysis. 


Cn  (MPa) 

C,  (MPa) 

C2  (MPa) 

ai  J 

a2 

136.5 

-135.2 

-1.270 

-45.28 

-1254 

Visualization  of  the  exponential  series  curve  fit  to  experimental  data  can  be  extremely 
helpful.  Use  the  Chart  Wizard  to  plot  both  experimental  and  predicted  shear  stress  values 
against  strain.  Choose  the  XY  (Scatter)  Plot  option.  Format  the  experimental  data  to  be 
displayed  as  symbols  and  the  predicted  data  as  a  continuous  line.  (See  Fig.  D.2). 


Nonlinear  regression  is  accomplished  by  the  ToolsJSolver  command  sequence.  Solver 
nonlinear  optimization  is  accomplished  via  the  method  of  Generalized  Reduced  Gradient. 
With  the  Solver  window  active  make  the  following  dialog  box  entries. 


1. 

2. 

3. 

4. 


5. 


Select  the  sum  of  the  difference  squared  cell  as  the  Set  Target  Cell. 

Choose  to  minimize  the  value  of  that  cell  by  setting  Equal  to:  Min. 

In  the  By  Changing  Cells:  enter  (or  highlight)  the  cells  containing  the  five  curve  fit 


the  Subject  to  Constraints:  box  select  Add  and  set  the  cell  containing  the 
3dulus  calculated  by  the  derivative  of  the  shear  stress  strain  curve  equal  to  the  cell 
ntaining  the  modulus  calculated  from  micromechanics. 

the  Subject  to  Constraints:  box  select  Add  and  set  the  predicted  stress  value 
responding  to  zero  strain  equal  to  zero. 
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6.  Select  the  Options  button  and  another  window  will  open. 

7.  In  the  Solver  Options  window  set  Tolerance:  to  0.1% 

8.  Make  sure  that  the  Use  Automatic  Scaling  option  is  turned  on 

9.  Click  OK  to  close  the  Solver  Options  window. 

1 0.  Click  Solve  to  begin  the  Solver  optimization. 

See  Fig.  D.2  for  an  sample  Solver  window  and  dialog  boxes. 


Fig.  D.2  Sample  Solver  window. 


Once  Solver  has  converged  on  a  solution,  it  is  recommended  using  the  new  curve  fit 
parameters  as  starting  values  for  another  Solver  iteration.  Continue  repeating  until  Solver 

no  longer  iterates. 


Better  curve  fits  of  MCT_77  analysis  results  against  experimental  data  with  shear  present 
can  be  achieved  if  the  incremental  modulus-strain  form  for  predicting  shear  stress  is  used 


=  (a_0  *  c_l  *  EXP(a_l  *  A1 1 )  +  a 


2  *  c_2  *  EXP(a_2  *  A1 1))*  (A12  -  A1 1)  +  C8 . 

(D.4) 
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Use  the  same  Solver  optimization  procedure  as  outlined  above.  Note  that  the  constant 
term  in  equation  (D.2)  does  not  appear  in  equation  (D.4).  An  additional  constraint  must 
be  added  so  that  equation  (D.l)  produces  zero  stress  given  zero  strain  when  optimizing 

equation  (D.4). 

Nonlinear  regression  can,  at  times,  require  some  trial  and  error  on  the  part  of  the  user  to 
arrive  at  a  suitable  solution.  If  the  trial  values  provided  in  Table  D.l  fail  to  produce  a 
convergent  solution,  use  the  graph  of  actual  and  predicted  values  as  a  guide  to  manually 
change  the  curve  fit  parameters  to  find  different  starting  values. 

D.2.  Curve  Fitting  Composite-Constituent  Relationships 

The  assumption  of  linear  elasticity  coupled  with  a  nonlinear  shear  stress-strain  curve 
implies  that  both  composite  and  constituent  material  properties  change  when  a  shear  load 
is  present.  MCT  computes  a  tangent  modulus  from  the  shear  stress-strain  curve  fit 
equation.  In  order  to  maintain  a  consistent  relationship  between  composite  and 
constituents,  the  constituents  must  reflect  the  changing  composite  shear  modulus: 
Changing  constituent  properties,  in  turn,  affect  other  composite  elastic  constants  besides 
shear  modulus.  A  set  of  second  order  polynomial  (quadratic)  curve  fits  defining  the 
composite-constituent  relationships  are  developed  in  a  one-time  operation  outside  the 
MCT  finite  element  code  using  micromechanics  and  an  Excel  spreadsheet. 

Although  it  is  the  composite  shear  strain  that  initiates  recalculation  of  all  material 
constants  within  the  code,  the  relationships  are  developed  by  varying  matrix  properties  in 
a  prescribed  manner  and  using  micromechanics  to  derive  new  composite  properties. 
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Once  a  set  of  constituent-composite  properties  are  defined,  the  relationships  are  defined 
in  terms  of  change  in  shear  modulus. 

To  characterize  the  matrix  for  micromechanics  analysis,  it  is  assumed  that  the  matrix  bulk 
modulus  is  independent  of  stress  state.  The  shear  modulus  is  varied  from  0%-100% 
(usually  in  increments  of  10%)  of  its  in  situ  value  as  determined  from  the  material  test 
suite  and  micromechanics.  Note  that  a  value  near  zero  (-0.1%)  is  used  instead  of  zero  to 
avoid  numerical  difficulties.  Having  defined  two  material  constants,  any  remaining 
matrix  properties  can  be  determined  using  standard  isotropic  material  relationships. 

The  results  of  each  discrete  increment  in  matrix  material  constants  provides  input  to  the 
micromechanics  model  for  calculating  corresponding  composite  constants.  The  matrix 
and  composite  material  constants  are  placed  in  an  Excel  spreadsheet.  Additionally,  four 
of  the  composite  material  constants  (E22,  G12,  G23,  V12)  and  the  matrix  shear  modulus  (G) 
must  be  expressed  as  a  ratio  of  their  current  value  to  the  initial  composite  shear  modulus. 
Examples  of  spreadsheet  layouts  that  effectively  present  the  above  analysis  can  be  found 
in  Tables  A.4,  B.5,  B.6,  C.6  and  C.7. 


The  graphing  routines  within  Excel  have  limited  curve  fitting  capability.  Use  the  Chart 
Wizard  to  create  four  plots  of  E22,  G,  G23,  and  v12  ratios  against  the  composite  shear 
modulus  ratio  for  each  increment.  Choose  the  XY  (Scatter)  Plot  option.  By  default,  the 
data  points  will  be  displayed  as  symbols.  Within  each  chart,  highlight  the  data  points 
(left  mouse  click  on  a  data  point)  and  then  right  mouse  click  for  a  menu.  Choose  Add 
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Trendline.  Under  the  Type  tab,  choose  Polynomial.  By  default  the  Order  will  be  set  to 
two.  Under  the  Options  tab,  check  the  dialog  boxes  for  Displaying  equation  on  chart  and 
Display  R-squared  value  on  chart.  Click  OK  to  close  the  Add  Trendline  window.  The 
quadratic  curve  fit  equation  for  each  of  the  four  plots  should  be  displayed.  Make  sure  the 
equation’s  ordinate  intercept  has  a  positive  or  zero  value.  If  it  is  negative,  highlight  the 
trendline  and  right  mouse  click  for  a  menu.  Choose  Format  Trendline.  Under  the 
Options  tab  check  the  dialog  box  for  Set  intercept  =.  The  default  value  is  zero.  Click  OK 
to  close  the  Format  Trendline  window. 


Coefficients  of  the  these  quadratic  equations  are  inputs  for  the  sdata  file  descnbed  in 
Appendix  F  Examples  of  material  ratio  plots  can  be  found  in  Figures  A.6-A.9,  B.7-B.14, 

and  C.7-C.14. 
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APPENDIX  E 

USING  EXCEL  FOR  OPTIMIZATION  OF  FAILURE  PARAMETERS 


The  MCT  failure  criteria  described  in  Section  4  requires  determination  of  five  failure 
parameters,  A,f,A4  f,  A2  m,  Aj m,  and  A4  m.  These  parameters  are  functions  of  ultimate 
constituent  strengths,  %{,  Sj,  %?,  %?,  SI2m,  and  Constituent  strengths  are 

derived  from  experimentally  determined  composite  lamina  ultimate  strengths. 
Determining  which  constituent  precipitates  composite  failure  is  necessary  for  establishing 
accurate  constituent  failure  values.  Identifying  which  constituent  precipitated  failure  in 
longitudinal  and  transverse  lamina  tension  and  compression  tests  is  intuitive  and 
straightforward.  Identifying  the  constituent  leading  to  shear  failure  is  more  problematic, 
especially  in  view  that  the  ultimate  constituent  shear  strengths  are  the  sole  contributors  to 
failure  parameters  As,  and  A4. 

This  appendix  outlines  a  procedure  for  optimizing  the  failure  parameters  by  using  non¬ 
linear  regression  analysis  in  a  manner  similar  to  that  used  in  Appendix  D  to  characterize 
shear  stress-strain  curves. 

E.l.  Failure  Criterion 

As  shown  in  Section  4,  the  general  form  of  the  stress  interactive  failure  criterion  used 
here  is 

Atl!  +  Ajl+A>Ii+A,l,=  1.  (E.D 
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Following  Hashin’s  delineation  of  two  failure  modes,  fiber  influenced  and  matrix 
influenced,  we  develop  a  form  of  (E.l)  for  each  constituent.  The  fiber  failure  criteria  has 

the  form 

±A'(l'J  +A{I{  =1,  (E-2) 

where  the  failure  parameters  were  determined  to  be 

1  a(  1 

±Af  = - ,  and  A  =7  ,  \2  • 

(*«)  fe) 

The  symbol,  ±,  indicates  that  the  appropriate  tensile  or  compressive  strength  value  be 
used  depending  on  the  stress  state. 


The  matrix  failure  criteria  has  the  form 


±a;(±i”)2  +a;i?  +a:i:  =\, 

(E.3) 

where  the  failure  parameters 

were  determined  to  be 

w  a 

II 

K> 

"ST  ^ 

a  3  ^ 

ro 

(E.4) 

1 

A  = 

4  .ter 

(E.5) 
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±  a  m 

i  i  (•sg-f+brri 

2teJ  J’ 
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1  {\  -  4m(±/mo)) 

=  ( 
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E.2.  Optimization  of  the  Failure  Parameters 

Microsoft®  Excel,  in  conjunction  with  the  Solver  analysis  option,  offers  a  convenient 
environment  for  implementing  a  nonlinear  regression  approach  to  solve  for  the  failure 
parameters.  Microsoft®  Excel  Solver  is  an  Add-in  analysis  package  included  with  the 
standard  spreadsheet  routines.  If  Solver  doesn’t  appear  under  the  Tools  menu  on  the 
toolbar,  then  it  must  be  installed  by  checking  the  Solver  option  box  under  Iools|Add- 

Ins|Available  Add-Ins. 

Data  from  off-angle,  balanced,  symmetric  laminates,  [±0]s,  provide  an  excellent  basis  for 
determining  optimized  failure  parameter  values.  These  laminates  produce  varying 
degrees  of  shear  and  normal  stresses  and  tend  to  fail  in  modes  allowing  analytical 
identification  of  the  constituent  precipitating  laminate  failure.  Other  off-angle  laminates, 
[+0]n  or  [-0]n  also  produce  varying  shear  and  normal  stresses  but  should  be  used  with 
caution  because  of  non-uniform  stress  states  caused  by  warping  of  the  test  specimen.  It  is 
important  to  generate  both  tension  and  compression  results  for  the  off-axis  laminates. 

From  ultimate  composite  failure  strengths  for  each  off-angle  laminate,  MCT_77  is  used 
to  calculate  the  constituents’  transversely  isotropic  stress  invariants  at  the  experimentally 
determined  point  of  composite  failure. 

Next,  failure  modes  for  each  laminate  tested  must  be  identified.  Composite  failure  is 
generally  fiber  dominated  at  lower  fiber  orientation  angles  (approximately  20°  and  less) 
and  matrix  dominated  at  higher  fiber  angles.  It  is  difficult  to  determine  failure  modes 
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visually,  so  it  becomes  more  convenient  to  do  so  analytically.  Good  initial  guesses  for 
constituent  ultimate  shear  strengths  can  be  derived,  via  MCT_77,  from  experimentally 
determined  ultimate  composite  shear  strengths  or  taken  from  the  literature  if  available. 
An  MCT  77  failure  analysis  of  each  off-axis  test,  using  the  initial  guesses,  tentatively 
identifies  the  failure  mode.  After  the  first  trial  optimization,  new  derived  values  for  As 
and  A4  be  can  used  to  check  the  failure  modes.  If  modes  change,  another  iteration  is 

required. 

Actual  optimization  begins  by  sorting  the  data  by  failure  mode  and  placing  the 
transversely  isotropic  stress  invariant  data  for  each  test  in  adjacent  columns  on  separate 
spreadsheets.  The  data  are  further  sorted  by  tension  or  compression  failures.  Cells  are 
designated  to  hold  the  failure  parameters  to  be  determined.  For  clarity,  labels  identifying 
the  parameters  are  placed  next  to  the  appropriate  cells.  It  is  also  very  advantageous  to 
use  the  Insert|Name|Defme  command  sequence  from  the  toolbar  to  name  each  cell  for  the 
appropriate  parameter.  Subsequent  equations  using  the  parameters  can  be  written  in  a 
more  natural  and  easier  to  read  manner.  See  Fig.  E.l  for  an  example  spreadsheet. 

The  value  of  A2,  as  determined  by  equation  (E.6),  depends  on  the  stress  invariants  at  the 
ultimate  transverse  strength,  I3m0  and  / 7°.  These  values  are  also  placed  on  the 
spreadsheet  and  named.  Note  that  A2  takes  on  different  tensile  and  compressive  values  as 

do  the  invariants  that  compose  it. 
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Fig.  E.l  Excel  spreadsheet  setup  for  determining  matrix  failure  parameters. 


i±jr 


41 


For  the  Solver  optimizer  to  work  properly,  all  invariants  must  be  scaled  to  a  reduced 
magnitude.  For  the  example  shown  above  using  SI  units,  the  invariants  are  scaled  by 
(MPa)2  or  1.0  x  1012.  The  failure  parameters  determined  by  the  optimization  will  be 

reduced  by  1.0  x  106. 


Next,  the  failure  state,  as  determined  by  the  failure  criterion  (in  this  example  the  matrix 
equation  (E.3)),  is  typed  into  the  first  cell  of  column  adjacent  to  the  invariant  data.  In  the 
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formula,  the  cell  reference  notation  ,  e.g.,  “B7’\  is  used  m  place  of  the  actual  invanant 
value  from  that  data  point.  For  example,  the  formula  might  look  like 

=  A_2t  *  B7  +  C7/S_23A2  +  D7/S_12A2 .  (E.7) 

The  formula  defined  for  the  first  cell  can  then  be  copied  and  pasted  to  the  remaining 
cells. 


The  next  adjacent  column  will  contain  the  square-of-the-differences  between  the 
predicted  failure  state  and  unity  (the  value  designating  failure).  Using  cell  reference 
notation,  the  formula  might  look  like 

=  (E7  -  F7)A2 .  (E-8) 

Again,  the  formula  can  be  written  for  the  first  cell  in  the  column  and  copied  and  pasted  to 
the  remaining  ones. 

At  the  bottom  of  the  column  containing  the  square-of-the-differences,  the 
Insert|Function|Math&T rig|Sum  command  sequence  is  used  to  sum  both  the  tension  and 
compression  cells  in  the  column.  This  cell,  designated  the  residual,  becomes  the  target  of 
subsequent  Solver  optimization  operations  to  minimize  it. 

It  appears  from  the  limited  analysis  done  within  this  research  that  the  more  data  points 
used  in  the  nonlinear  optimization,  the  better  the  convergence  to  optimized  failure 
parameters.  Convergence  is  further  enhanced  by  constraining  optimization  to  achieve 


E-6 


NS  WCCD-65-TR- 1 999/1 5 


failure  states  of  unity  for  the  limiting  cases  of  pure  transverse  tension  and  compression 
loads,  and  positive  failure  states  for  all  other  cases. 


Nonlinear  regression  is  accomplished  executing  the  loolsjSolver  command  sequence. 
Solver  nonlinear  optimization  is  accomplished  via  the  method  of  Generalized  Reduced 
Gradient.  With  the  Solver  window  active  make  the  following  dialog  box  entries: 

1 .  Select  the  sum  of  the  difference  squared  cell  as  the  Set  Target  Cell. 

2  Choose  to  minimize  the  value  of  that  cell  by  setting  Equal  to:  Min. 

”3.  In  the  By  Changing  Cells:  enter  (or  highlight)  the  cells  containing  the  five  curve 

fit  parameters.  .  .  .  . .. 

4.  In  the  Subject  to  Constraints:  box  select  Add.  Set  the  cell  contaimng  the  failure 
states  for  the  transverse  only  loads  to  unity  and  all  other  failure  states  to  greater 
than  or  equal  to  zero  (positive  definite).  . 

5.  In  the  Subject  to  Constraints:  box  select  Add  and  set  the  predicted  stress  value 
corresponding  to  zero  strain  equal  to  zero. 

6.  Select  the  Options  button  and  another  window  will  open. 

7.  In  the  Solver  Options  window  set  Tolerance:  to  0. 1  % 

8.  Make  sure  that  the  Use  Automatic  Scaling  option  is  turned  on 

9.  Click  OK  to  close  the  Solver  Options  window. 

1 0.  Click  Solve  to  begin  the  Solver  optimization. 


See  Fig.  E.2  for  an  sample  Solver  window  and  dialog  boxes. 


The  results  of  a  nonlinear  analysis  regression  are  highly  dependent  on  starting  values.  As 
mentioned  previously,  ultimate  constituent  shear  strength  values  derived  from 
experimentally  determined  composite  results  provide  excellent  starting  values.  There 
appear  to  be  numerous  local  minima  in  the  failure  surface  defined  by  all  possible 
(realistic  or  not)  failure  parameters.  Hence,  it  is  important  to  do  numerous  optimizations 
trials  using  different  starting  points  to  determine  the  global  minima.  Rows  21  through  25 
in  the  example  spreadsheet  show  the  results  of  five  trials. 
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Fig.  E.2  Sample  Solver  window. 


Bad  data  points  will  skew  the  nonlinear  regression  procedure  and  produce  less  than 
optimum  failure  parameters.  Bad  data  points  can  be  identified  by  their  failure  to 
converge  to  near  unity  as  shown  in  row  8  of  Fig.  E.l.  In  this  case,  simply  exclude  the 
contribution  of  that  data  point’s  square-of-the-difference  to  the  residual  calculation. 

Once  the  global  minimum  has  been  found,  the  values  for  the  ultimate  constituent  shear 
strengths  are  placed  in  the  fdata  file  for  MCT_77  failure  analysis.  Note  that  the  example 
shown  above  was  an  optimization  for  matrix  parameters  only.  A  similar  procedure  must 


be  used  to  determine  A4  for  the  fiber. 
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APPENDIX  F 

MCTJ77  USER’S  GUIDE 


F.l.  Introduction 

MCT  77  is  a  traditional  displacement-based,  linear  elastic,  finite  element  code.  What 
makes  MCT_77  different  from  other  finite  element  codes  is  the  implementation  of 
Multicontinuum  Theory  (MCT).  MCT  allows  the  user  to  extract  constituent  stress  and 
strain  fields  from  those  of  the  composite  during  the  course  of  a  routine  structural  analysis. 

The  MCT  77  finite  element  software  code  is  written  in  FORTRAN  77  for  maximum 
portability.  It  is  intended  as  a  research  code  and  is  not  optimized  for  speed  as  a 
commercial  program  might  be.  The  code  is  extensively  commented  and  written  in  a 
modular  fashion  for  maximum  readability  and  to  facilitate  changes.  It  contains 
approximately  seventy  subroutines  and  ten  thousand  lines  of  code.  It  has  been  compiled 
and  run  under  DEC,  IBM,  and  Microsoft  versions  of  FORTRAN.  The  code  has  been 
tested  on  a  DEC  ALPHA  workstation,  an  IBM  workstation,  and  a  Gateway  Personal 

Computer. 

The  MVT  77  program  was  influenced  by  the  following  books  and  their  example  finite 
element  codes  contained  within  references  [1, 2, 3] 

In  particular,  Reference  [1]  was  used  as  the  basis  for  many  of  the  finite  element  related 
formulations  encoded  as  various  subroutines  in  MCT_77. 
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A  translator,  ANS2MCT,  has  been  written  allowing  use  of  the  pre-processing  capabilities 
of  ANSYS  to  create  the  geometry  and  finite  element  mesh  of  the  structure  being 
analyzed.  If  requested,  MCT_77  writes  results  files  that  can  be  read  by  ANSYS  for  post¬ 
processing. 

An  MCT_77  analysis  is  defined  by  three,  four,  or  five  user-generated  ASCII  format, 
input  files:  fdata  (failure  data),  gdata  (grid  data),  ldata  (load  data),  mdata  (material 
data),  and  sdata  (nonlinear  shear  data).  Fdata  and  sdata  files  are  required  only  if  the  user 
requests  a  failure  analysis  and/or  wants  to  include  nonlinear  shear  effects. 

Data  in  each  input  file  is  read  sequentially.  Order  of  the  input  lines  is  fixed  and  all  input 
files  share  a  common  format.  Each  line  of  input  is  preceded  by  a  line  of  descriptive  text. 
While  some  lines  of  input  may  depend  on  the  status  of  flags  set  by  the  user,  the  lines  of 
descriptive  text  must  always  be  included.  Within  the  User’s  Guide,  each  line  of  input  is 
referenced  to  a  corresponding  line  of  descriptive  text  which  is  designated  by  the  LINE 
notation  and  a  sequence  number. 

Generally,  it  is  easier  to  use  the  ANSYS  pre-processor  to  create  new  gdata  and  ldata  files 
and  modify  existing  mdata,  fdata,  or  sdata  files  to  fit  a  new  case  rather  than  create  new 
ones  from  scratch. 

The  user  should  be  aware  that  the  MCT  77  software  contains  virtually  no 
capability  to  detect  or  correct  contradictory  or  erroneous  data  input. 
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The  analysis  results  are  printed  in  the  odata  (output  data)  file.  The  odata  file  echos  all 
input  data  and  prints  analysis  results  requested  by  the  user  in  the  ldata  file. 


The  remainder  of  this  appendix  describes  each  of  the  five  input  files  and  the  ANSYS-to- 
MCTJ77  translator  in  detail.  Each  input  file  description  uses  the  following  framework: 


1  LINE  #  of  the  descriptive  text  delineating  a  section  of  the  file, 

2.  FORMAT  statement  (FORTRAN  77  convention)  delineating  the  alphanumeric 
FIELDS  defining  analysis  variables  that  form  the  LINE, 

3.  FIELD  #  provides  a  textual  definition  of  each  variable, 

4.  OPTIONS  available  for  each  FIELD,  and 

5.  NOTES  containing  supplementary  information  for  the  LINE  or  FIELD. 
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F.2.  ANSYS-to-MCT_77  Translator 

A  translator,  ANS2MCT,  has  been  written  allowing  use  of  the  pre-processing  capabilities 
of  ANSYS  to  create  the  geometry  and  finite  element  mesh  of  the  structure  being 
analyzed.  The  translator  converts  ANSYS  database  information  to  MCT_77  readable 
format.  ANS2MCT  translates  the  following  elements: 


Table  F.l  ANSYS  to  MCT  element  translation. 


ANSYS  Element  Type 

MCT  Element  Type _ 

PLANE42 

QUAD4 

PLANE82 

OUAD8 

SOLED45 

BRICK8 

SOLID46 

BRICK8 

SHELL93 

SHELL8 

SHELL99 

SHELL8 

Note  that  not  every  capability  of  the  ANSYS  elements  listed  in  Table  F.l  is  translated  to 
MCTJ77.  The  user  should  manually  verify  each  translation.  ANS2MCT  translates  the 

following  element  attributes: 

1 .  Nodal  coordinates, 

2.  Nodal  connectivity, 

3.  Plane  stress/strain  option, 

4.  Material  number, 

5.  Non-conforming  element  formulation  option,  ; 

6.  Local  coordinate  system  (ANSYS)  <e>  Principal  material  direction  (MCT_77), 

7.  Number  of  layers,  and 

8.  Layer  thickness. 


ANS2MCT  translates  the  following  analysis  options: 

1 .  Number  of  load  steps  (always  set  to  one), 

2.  Analysis  type  (always  set  to  linear  elastic),  and 

3.  Failure  analysis  (always  set  to  OFF). 
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ANS2MCT  defaults  to  the  following  output  options: 

1 .  No  printing  of  field  results  in  odata, 

2.  No  plot  results, 

3  Stress/strain  values  are  interpolated  to  the  nodes, 

4.  Stress/strain  values  are  referenced  to  the  principal  material  coordinate  system 
(Required  for  ANSYS  post-processing  compatibility),  and 

5.  Write  ANSYS  readable  results  files  is  ON. 


ANS2MCT  translates  the  following  boundary  conditions: 

1 .  Prescribed  homogeneous  nodal  displacements, 

2.  Prescribed  non-homogeneous  nodal  displacements, 

3.  Nodal  forces, 

4.  Uniform  temperature, 

5.  Constraint  equations,  and 

6.  Element  pressures. 

MCTJ77  writes  ANSYS  readable  results  files  for  graphical  post-processing.  The  results 
files  can  be  written  for  single  or  multi-step,  linear-elastic  or  nonlinear-elastic,  failure  or 
non-failure,  MCT_77  analyses.  The  following  results  files  are  created: 

1 .  mct2ans.c  (composite  resuits), 

2.  mct2ans.r  (reinforcement  results), 

3.  mct2ans.m  (matrix  results),  and 

4.  damage.dat  (composite  damage  state  if  failure  analysis). 

The  mct2ans  files  for  the  composite  contain  nodal  displacements,  stresses,  and  strains  for 
each  load  step.  The  mct2ans  files  for  each  constituent  contains  stress  and  strain 
information  for  each  load  step.  The  damage.dat  file  contains  material  numbers 
corresponding  to  the  composite  damage  state  at  each  gauss  point  in  the  model. 

The  MCT_77  finite  element  program  and  the  ANS2MCT  translator  can  be  initiated 
within  the  ANSYS  graphical  user  interface  (GUI)-  Five  ANSYS  macros  that  execute 
MCT  77  related  functions  have  been  incorporated  as  pushbuttons  on  the  ANSYS  toolbar. 


F-5 


NSWCCD-65-TR- 1999/1 5 


l,A>  AMSYS  Tnnlbar 
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Fig.  F.l  ANSYS  toolbar  with  MCTJ77  translator. 


The  pertinent  MCT  toolbar  buttons  are: 


1. 

2. 

3. 

4. 

5. 

6. 


WRITEMCT:  Create  MCT_77  input  files  that  defining  the  finite  element  mesh 
(gdata)  and  the  boundary  conditions  (ldata), 

MCT:  Executes  the  MCT_77  finite  element  solver, 

LAMINA:  Reads  MCT_77  composite  results  a  specified  lamina  for  all  elements, 
FIBER:  Reads  MCT_77  fiber  results  for  a  specified  lamina  for  all  elements, 
MATRIX:  Reads  MCT_77  matrix  results  for  a  specified  lamina  for  all  elements,  and 
DAMAGE:  Reads  MCT_77  material  numbers  for  a  specified  lamina  for  all  elements. 


WRITEMCT  macro  writes  the  ANSYS  database  to  a  text  file  and  executes  ANS2MCT 
which  translates  the  database  information  to  gdata,  and  ldata  files.  The  user  is 
responsible  for  any  additional  files  required  for  a  MCT_77  analysis.  At  a  minimum,  the 
user  must  generate  a  mdata  file  to  define  material  properties  for  the  finite  element  model. 
Note  that  the  material  number  assigned  by  ANSYS  to  each  element  must  match  the 
material  number  assigned  by  the  user  in  the  mdata  file. 


MCT  macro  runs  the  MCT_77  solver  based  on  input  files  created  by  ANS2MCT  and  the 
user.  WRITEMCT  and  MCT  macros  are  not  combined  giving  the  user  an  opportunity  to 
check  ldata  and  gdata  files  and  create  an  associated  mdata  file. 
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Currently,  an  ANSYS  solution  must  be  executed  for  the  model  before  storage  positions 
are  created  in  the  ANSYS  database  for  displaying  MCT_77  post-processing  results. 
LAMINA,  FIBER,  MATRIX,  and  DAMAGE  macros  read  MCT_77  results  files  into  the 
ANSYS  database.  The  user  will  be  queried  for  a  specific  load  step  (defaults  to  load  step 
1).  ANSYS  does  not  have  the  capability  to  store  stress  and  strain  values  for  the 
composite  and  both  constituents.  Each  time  a  macro  reads  an  MCT_77  results  file,  it 
overwrites  any  existing  displacement,  stress,  and  strain  values.  Due  to  MCT-ANSYS 
interface  difficulties,  lamina  results  are  read  in  one  at  a  time.  The  results  are  read  into 
layer  0  (zero)  position  in  the  ANSYS  database.  Subsequent  reading  of  MCT_77  results 
overwrites  data  stored  in  this  location  also.  DAMAGE  data  are  read  into  an  element  table 
for  listing  or  plotting. 
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F.3.  Fdata  File 

The  fdata  file  contains  ultimate  strength  data  for  the  composite  and  its  constituents  (fiber, 
and  matrix).  Currently,  failure  in  only  one  composite  material  can  be  modeled.  Fdata 
and  mdata  files  must  be  structured  in  a  prescribed  manner  to  conduct  a  failure  analysis. 

The  mdata  file  must  be  have  the  following  format: 

1.  Material  #1  data  must  contain  undamaged  composite  material  constants. 

2  Material  #2  data  must  contain  material  #1  constants  assuming  damaged  matrix. 

3‘  Material  #3  data  must  contain  material  #1  constants  assuming  damaged  fiber. 

4.  Materials  #4-#7  data  contain  constituent  material  properties  corresponding  to 
Materials  #l-#3,  i.e.,  undamaged  and  damaged  values  for  the  fiber  and  matrix. 

The  fdata  file  must  have  the  following  format: 

1 .  Material  #1  data  must  correspond  to  material  #1  in  the  mdata  file. 

2.  Material  #2  data  must  correspond  to  the  undamaged  fiber  material  in  mdata. 

3.  Material  #3  data  must  correspond  to  the  undamaged  matrix  material  in  mdata. 

Changing  material  properties,  inherent  in  failure  analyses,  requires  incremental 
application  of  the  load  and  iterations  within  each  load  step  to  achieve  equilibrium.  Two 
tolerances,  deflection  (DTOL)  and  load  (FTOL)  control  equilibrium  iterations.  The  user 
should  monitor  convergence  of  failure  analyses  carefully  and  modify  tolerances  (default 
values  are  DTOL=7.0  and  FTOL=1.0E-4)  as  required.  A  sample  fdata  file  is  shown  in 
Fig.  F.2.  The  file  contains  six  required  LINES.  A  LINE-by-LINE  description  follows. 
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‘Failure  data  title 
I  Boron/5505,  50%  FVF 
[Convergence  tolerances,  deflection,  load 
1.0D0  1.0D-3 

‘Material  Sets 

'Material, Composite, SI  1T,S22T,S33T,S12T,S13T,S23T,(repeat  for  compress) 

1379.0D6  86.6D6  00.0D6  125.5D6  00.0D6  00.0D6 

-2861. 0D6  -262.0D6  00.0D6  125.5D6  00.0D6  00.0D6 

‘Material, Reinforcement, SI  1T.S22T ,S33T ,S1 2T.S1 3T.S23T, (repeat  for  compress) 
2708.0D6  00.0D6  00.0D6  166.5D6  00.0D6  00.0D6 

-5620  0D6  00.0D6  00.0D6  166.5D6  00.0D6  00.0D6 

‘Material, Matrix, S1 1  T,S22T,S33T,S1 2T.S1 3T.S23T, (repeat  for  compress) 

138.0D6  69.8D6  6.68D6  84.5D6  25.9D6  28.2D6 

-345.0D6  -208.8D6  -20.0D6  84.5D6  -77.4D6 _ 28.2D6_ - 


Fig.  F.2  Sample  /data  file. 


LINE  1.  *Failure  data  title 
Number  of  fields:  1 

FORMAT:  (A60)  #  „ 

FIELD  1:  General  description  of  the  failure  data  defined  in  the  fdata  tile. 

OPTIONS:  None  .  . 

NOTE:  Descriptor  is  printed  as  an  identifier  at  the  beginning  of  the  odata  tile. 

LINE  2.  *  Convergence  tolerances,  deflection,  load. 

Number  of  fields:  2 

FORMAT:  (IX, 6G12.7)  „ ,  ^  ,  tx2 

FIELD  1 :  Maximum  change  allowed,  within  a  load  step,  of  the  displacement  L  norm 

before  failure  of  the  structure  is  declared  and  program  execution  ends. 

OPTIONS:  None. 

NOTE:  Displacement  tolerance  is  defined  as: 

Delta  displacement^  for  an  equilibrium  iteration 
Displacement  for  the  first  equilibrium  iteration 

It  is  easiest  to  think  of  this  tolerance  as  a  percentage  change  in  displacement.  Except  for 
the  first  iteration  in  the  load  step,  the  percentage  change  in  displacement  norm  is  much 
less  than  100%  and  should  decrease  as  equilibrium  is  achieved.  When  failure  occurs 
during  equilibrium  iteration,  the  structure  can  no  longer  support  the  load  and 
displacements  grow  without  bound.  The  displacement  tolerance  stops  the  program  before 
numeric  singularities  occur.  Typically,  values  around  6.0  work  well  (7.0  is  default). 
FIELD  2:  Specifies  the  maximum  load  tolerance  defined  as  a  ratio  of  the  L  norm  of  the 
internal  load  vector  norm  to  the  L2  norm  of  the  external  load  vector  before  equilibrium  is 

satisfied. 

OPTIONS:  None. 
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NOTE:  Changing  material  properties,  due  to  failure  or  nonlinear  shear,  during  an 
incremental  analysis  cause  the  external  and  internal  (resisting)  load  vectors  to  be  out  o 
balance  The  load  tolerance  controls  the  number  of  times  that  the  delta  load  difference 
between  the  two  load  vectors  is  applied  to  the  structure  before  equilibrium  is  achieved 

(1.0E-4  is  the  default). 


LINE  3.  *Material  Sets 
Number  of  fields:  1 

FORMAT:  (IX,  915)  .  , 

FIELD  1 :  Number  of  materials  for  which  failure  data  is  to  be  read. 

OPTIONS:  3  (min)  to  8  (max).  See  Note.  .  . 

NOTE:  The  code  uses  failure  data  for  only  one  composite  material  (corresponding  to  the 
first  material  listed  in  mdata)  and  its  associated  reinforcement  and  matnx  constituents. 
As  such,  only  a  value  designating  three  materials  (composite  and  two  constituents)  will 

produce  reliable  results. 


LINE  4.  *Material,Composite,S  1 1  T,S22T,S33T,S  1 2T,S  1 3T,S23T, (repeat  for  compress) 
Number  of  fields:  12 


FORMAT:  (IX, 6G12.7) 

NOTE  1 :  Two  lines  are  required,  1st  -  tension,  2  -  compression 
NOTE  2:  Composite  values  are  used  to  calculate  Tsai-Wu  failure  criterion. 

FIELD  1:  Ultimate  composite  tensile  strength  in  the  fiber  (+Sn)  direction. 

OPTIONS:  None  .  ....  ^  v- 

FIELD  2:  Ultimate  composite  tensile  strength  transverse  to  the  fiber  direction  (  022) m 

the  plane  of  the  lamina. 


OPTIONS:  None 

FIELD  3:  Ultimate  composite  tensile  strength  transverse 
of  the  plane  of  the  lamina. 


to  the  fiber  direction  CS33)  out 


OPTIONS:  None 

NOTE:  Currently  not  used.  . 

FIELD  4:  Ultimate  composite  shear  strength  in  the  12  direction  (Si 2). 

OPTIONS:  None  .  . 

FIELD  5:  Ultimate  composite  shear  strength  in  the  13  direction  (Si 3). 


OPTIONS:  None 
NOTE:  Currently  not  used. 

FIELD  6:  Ultimate  composite  shear  strength  in  the  23  direction  (S23T 
OPTIONS:  None 


FIELDS  7-12:  Repeats  FIELDS  1-6  for  compression  values. 

LINE  5.  *Material,Reinforcement,Sl  lT,S22T,S33T,S12T,S13T,S23T,(repeat  for 
compress) 

Number  of  fields:  12 
FORMAT:  (1X,6G12.7) 

NOTE:  Two  lines  are  required,  1st-  tension,  2  -compression  + 

FIELD  1:  Ultimate  reinforcement  tensile  strength  along  the  fiber  axis  (  Sn). 
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OPTIONS:  None 

FIELD  2:  Ultimate  reinforcement  tensile  strength  transverse 


to  the  fiber  axis  CS22)  in  the 


plane  of  the  lamina. 

OPTIONS:  None 

NOTE:  Currently  not  used.  .  f 

FIELD  3:  Ultimate  reinforcement  tensile  strength  transverse  to  the  fiber  axis  (  S33) 

the  plane  of  the  lamina. 

OPTIONS:  None 

NOTE:  Currently  not  used.  . 

FIELD  4:  Ultimate  reinforcement  shear  strength  in  the  12  direction  (Si 2). 


OPTIONS:  None  .  ,  .  t  ,c  , 

FIELD  5:  Ultimate  reinforcement  shear  strength  in  the  13  direction  (bo) 


OPTIONS:  None 

NOTE:  Currently  not  used.  . 

FIELD  6:  Ultimate  reinforcement  shear  strength  in  the  23  direction  (S23) 


OPTIONS:  None 
NOTE:  Currently  not  used. 

FIELDS  7-12:  Repeats  FIELDS  1-6  for  compression  values. 


LINE  6.  ‘Material, Matrix, si  lT,S22T,S33T,S12T,S13T,S23T,(repeat  for  compress) 

Number  of  fields:  12 
FORMAT:  (1X.6G12.7) 

NOTE:  Two  lines  are  required,  1st  -  tension,  2n  -  compression. 

FIELD  1:  Matrix  in-situ  stress  value  in  the  1 1  direction  (on)  when  the  lamina  is  loaded 
in  the  22  direction  to  its  ultimate  value  (822)- 
OPTIONS:  None 

NOTE:  Currently  not  used.  .  .+_  ,  .  , 

FIELD  2:  Ultimate  matrix  in-situ  tensile  strength  transverse  to  the  fiber  axis  (  822)  in  t  e 

plane  of  the  lamina. 

OPTIONS:  None 

NOTE:  This  is  backed  out  of  transversely  loaded  (S22)  lamina  test. 

FIELD  3:  Matrix  in-situ  stress  value  in  the  33  direction  (033)  when  the  lamina  is  loaded 
in  the  22  direction  to  its  ultimate  value  (S22)- 
OPTIONS:  None 

NOTE:  This  is  backed  out  of  transversely  loaded  (S22)  lamina  test. 

FIELD  4:  Ultimate  matrix  in-situ  shear  strength  in  the  12  direction  (Si 2). 

OPTIONS:  None 

NOTE:  This  is  backed  out  ofa  shear  loaded  (S12)  lamina  test. 

FIELD  5:  Ultimate  matrix  in-situ  tensile  strength  along  the  fiber  axis  (Si  1). 

OPTIONS:  None 

NOTE:  Currently  not  used.  . 

FIELD  6:  Ultimate  matrix  in-situ  shear  strength  in  the  23  direction  (S23)- 

OPTIONS:  None 

NOTE:  This  is  backed  out  of  a  shear  loaded  (S23)  lamina  test. 

FIELDS  7-12:  Repeats  FIELDS  1-6  for  compression  values. 
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F.4.  Gdata  File 

The  gdata  file  contains  infoimation  defining  elements  that  discretizes  the  structure  being 
analyzed.  A  sample  gdata  file  is  shown  in  Fig.  F.3.  The  gdata  file  contains  six  required 
LINES.  A  LINE-by-LINE  description  follows. 


*GRID  MODEL  TITLE 

SINGLE  ELEMENT,  8-NODE  BRICK  MODEL 
*GRID,  elem  type,  plain  stress/strain,  NCEF,  #elms,  #pts 
BRICK8  0  0  1  8 

‘LAYER,  #layer/elm,  thk 


4  1.D0 

COORDINATES,  define  global  nodal  coordinates. 

1  0.0000 

0.0000 

0.0000 

2  1.0000 

0.0000 

0.0000 

3  1.0000 

1.0000 

0.0000 

4  0.0000 

1.0000 

0.0000 

5  0.0000 

0.0000 

1.0000 

6  1 .0000 

0.0000 

1.0000 

7  1.0000 

1.0000 

1.0000 

8  0.0000 

1.0000 

1.0000 

‘CONNECTIVITY,  define  element  connectivity. 
112345678 

‘ORIENTATION,  define  element  material,  layer  principal  material  direction. 
1  1  30. -30. -30.  30. - - - - 

Fig.  F.3  Sample  gdata  file. 


LINE  1.  *GRID  MODEL  TITLE 
Number  of  fields:  1 
FORMAT:  (A60) 

FIELD  1  :General  description  of  the  mesh  defined  in  the  gdata  file. 
OPTIONS:  None 

NOTE:  Descriptor  is  printed  as  an  identifier  at  the  beginning  of  the  odata  file. 

LINE  2.  *GRID,  elem  type,  plain  stress/strain,  NCEF,  #elms,  #pts 
Number  of  fields:  5 
FORMAT:  (1X,A1 0,415) 
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FIELD  1:  Type  of  element  used. 

OPTIONS: 

a)  QUAD4  (4  node  quadrilateral), 

b)  QUAD8  (8  node  quadrilateral), 

c)  BRICK8  (8  node  brick), 

d)  SHELL8  (8  node  shell).  . 

NOTE:  Currently  only  one  type  of  element  can  be  use  within  an  analysis. 

FIELD  2:  Element  kinematic  assumptions. 

OPTIONS: 

a)  QUAD4  -  (1)  =  plane  stress  or  (2)  =  plane  strain, 

b)  QUAD8  -  (1)  =  plane  stress  or  (2)  =  plane  strain, 

c)  BRICK8  -  (0)  only  (folly  three  dimensional), 

d)  SHELL8-(1)  only  (plane  stress). 

FIELD  3:  Request  that  a  non-conforming  element  formulation  be  used. 

OPTIONS: 

a)  QUAD4  -  Off  (0)  or  On  (1), 

b)  QUAD8-  Off  (0)  only, 

c)  BRICK8  -  Off  (0)  or  On  (1), 

d)  SHELL8- Off  (0)  only. 

NOTE:  This  option  is  only  applicable  to  linear  elements. 

FIELD  4:  Number  of  elements  in  the  model. 

OPTIONS:  None 

FIELD  5:  Number  of  nodes  in  the  model. 

OPTIONS:  None 

LINE  3.  *LAYER,  #layer/elm,  thk 
Number  of  fields:  2 
FORMAT:  (1X,I5,G10.5) 

FIELD  1 :  Number  of  layers  per  element, 

OPTIONS:  1  (min)  to  26  (max)  layers. 

FIELD  2:  Element  thickness  in  the  C,  direction  as  determined  by  element  connectivity. 
OPTIONS: 

a)  QUAD4  -  constant  thickness, 

b)  QUAD8  -  constant  thickness, 

c)  BRICK8  -  constant  thickness.  Note  that  this  value  must  agree  with  the  thickness  as 
determined  by  distance  between  the  nodes  in  C,  direction. 

d)  SHELL8  -  constant  thickness. 

LINE  4.  ““COORDINATES,  define  global  nodal  coordinates. 

Number  of  fields:  Variable  3  to  4  depending  on  element  type. 

FORMAT:  (1X,I5,10G15.8) 

FIELD  1:  Node’s  numeric  label. 

OPTIONS:  None 

NOTE:  Nodes  must  be  numbered  sequentially  starting  with  1  having  no  gaps  or 
duplicates. 

FIELD  2:  Value  of  x  coordinate  in  global  Cartesian  coordinate  system. 
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OPTIONS:  None  ..  t 

FIELD  3:  Value  of  y  coordinate  in  global  Cartesian  coordinate  system. 

OPTIONS:  None  .  ..  t 

FIELD  4:  Value  of  z  coordinate  in  global  Cartesian  coordinate  system. 


OPTIONS:  ..  t  _ 

a)  QUAD4  -  elements  must  be  in  the  x-y  plane,  z  coordinate  set  to  zero, 

b)  QUAD8  -  elements  must  be  in  the  x-y  plane,  z  coordinate  set  to  zero, 

c)  BRICK8  -  three  dimensional  nodal  coordinates, 

d)  SHELL8  -  three  dimensional  nodal  coordinates. 

NOTE:  Number  of  lines  must  equal  number  of  nodes  specified  in  LINE  2. 


LINE  5.  CONNECTIVITY,  define  element  connectivity. 

Number  of  fields:  Variable  4-8  depending  on  element  type. 

FORMAT:  (IX, 915) 

FIELD  1:  Node’s  numeric  label. 

OPTIONS:  None  . .  . 

NOTE:  Nodes  must  be  numbered  sequentially  starting  with  1  having  no  gaps 

FIELDS  2-9:  Element  connectivity.  Order  is  counterclockwise,  comer  nodes  first, 
midside  nodes  next  (if  any). 

OPTIONS:  None  T  ..  TmIP, 

NOTE:  Number  of  lines  must  equal  number  of  elements  specified  in  LLNL  i. 


LINE  6.  *ORIENTATION,  define  element  material,  layer  principal  material  direction. 
Number  of  fields:  Variable  3-16  depending  on  number  of  layers  in  element. 

FORMAT:  (1X,2I5,14G5.1) 

FIELD  1:  Element  number. 

OPTIONS:  None  . 

NOTE:  Elements  must  be  numbered  sequentially  starting  with  1  having  no  gaps  o 


duplicates. 

FIELD  2:  Material  type  for  this  element. 

OPTIONS:  None  .  ,  ,  c , 

NOTE:  Only  one  material  type  per  element  is  allowed  regardless  of  the  number  o  ayer  . 

Must  correspond  to  one  of  eight  possible  material  types  defined  in  mdata. 

FIELDS  3-16:  Principal  material  direction  in  each  layer  as  an  angle  of  rotation  about 
positive  z  axis  (positive  counterclockwise)  of  the  element  coordinate  system.  See  No  es 

2  and  3. 

OPTIONS:  None  „  ,  .  _  , .  T  Txrc  o 

NOTE  1 :  Number  of  lines  must  equal  number  of  elements  specified  in  LIME.  2. 

NOTE  2:  Number  of  fields  must  equal  number  of  layers  specified  in  LINE  3. 

NOTE  3’  Each  element  has  an  associated  coordinate  system.  This  system  is  the  frame  ot 
reference  for  determining  principal  (fiber)  material  direction.  The  element  coordinate 
system  is  determined  by  a  projection  of  the  global  coordinate  system  (referred  to  as  the 
“shadow”  axis)  onto  the  element  surface  using  the  following  rules: 
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a) 

b) 

c) 

d) 

e) 

f) 


Element  normal  (z’  direction)  is  computed  using  the  cross  product  of  the  positive  £, 
and  11  vectors  originating  at  the  first  node.  Note  that  these  vectors  are  generally  not 

If  the  normal  ha^an  absolute  z  component  (global  reference)  greater  than  1  degree, 
the  shadow  coordinate  system  coincides  with  the  global  coordinate  system,  i.e.,  no 


rotation  is  necessary  for  projection.  uncmvPv 

If  the  normal  is  within  one  degree  of  the  global  x-y  plane  and  has  a  POSITIVE  y 
component,  rotate  the  shadow  coordinate  system  -90  degrees  about  the  global  x- 
axis  and  project  (positive  shadow  z  axis  coincides  with  positive  global  yaxis). 

If  the  normal  is  within  one  degree  of  the  x-y  plane  and  has  a  NEG AT  y 

component,  rotate  shadow  coordinate  system  +90  degrees  about  the  global  x-axis 
and  project  (positive  shadow  z  axis  coincides  with  negative  global  y  axis) 

If  the  normal  is  within  one  degree  of  the  global  x-y  plane  and  within  one  degree  of 
the  POSITIVE  x-axis,  rotate  the  shadow  coordinate  system  +90  degrees  abou  e 
global  y-axis  and  project  (positive  shadow  z  axis  coincides  with  positive  global  x 


If  the  normal  is  within  one  degree  of  the  global  x-y  plane  and  within  one  degree  of 
the  NEGATIVE  x-axis,  rotate  the  shadow  coordinate  system  -90  degrees  about  the 
global  y-axis  and  project  (positive  shadow  z  axis  coincides  with  the  negative  global 


x  axis). 
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F.5.  Ldata  File 

The  ldata  file  contains  information  defining  boundary  conditions  of  the  structure  being 
analyzed.  A  sample  ldata  file  is  shown  in  Fig.  F.4.  The  ldata  file  contains  thirteen 
required  LINES.  A  LINE-by-LINE  description  follows. 


•LOAD  CASE  TITLE 
Combined  boundary  conditions,  single  brtck8  element. 

•LOAD  CONTROL,  #  steps  total,  #steps  to  execute,  analysis  type 

•LOAD^NDEX,  #BCs,  #N-forces,  #N-displ,  #MPCs,  #E-Press,  B-force,  T-loads 

7  4  14  10  1 

•FIXED  BC,  define  homogenous  nodal  dot  displacements. 

1  111 
2  011 

3  001 

4  101 

5  110 

6  010 

8  100 

•NODAL  FORCES,  define  nodal  dof  forces. 

5  00.00D0  00.00D0  100.00D6 

6  00.00D0  00.00D0  100.00D6 

7  OO.OODO  00.00D0  100.00D6 

.  8  OO.OODO  OO.OODO  100.00D6 

rNODAL  DISPS,  define  nonzero  nodal  dof  displacements. 

3  021  0.0000  1.0D-3  0.0000 

|*MPC,  define  multi-point  dof  contraints. 

2  0.00000E+00  3  2  1.0000  4  2  -1.0000 

2  0  00000E+00  3  2  1.0000  7  2  -1.0000 

2  O.OOOOOE+OO  3  2  1.0000  8  2  -1.0000 

•ELM  PRESSURES,  define  element  edge  or  face  pressures. 

1  3  -100. 0D6  -100. 0D6  -100.0D6  -100.0D6 

rBODY  FORCE,  define  forces  due  to  accelerations. 
pTEMPERATURE  LOAD,  define  uniform  or  nodal  temperatures 

•OUTPUT  CONTROL,Maxprn,prnsteps,Local,Gauss,Lay,Fail,ANSYS,Plot 
•PRINT  CONTROL, Inc, Displs,Nsts,Nstn,Asts,Astn,Ests,Estn,Psts,lsts,Rfor 

1 YNNYYNNNYN 

•PLOT  CONTROL, {force  id, force#, blank, blank}  or  {field  id,node,lay,elem} 

23  1  0  0  3  2  1  1 _ _ _ _ _ 


Fig.  F.4  Sample  ldata  file. 
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LINE  1.  *LOAD  CASE  TITLE 
Number  of  fields:  1 


FORMAT:  (A60)  .  .  ,  . 

FIELD  1  :General  description  of  the  boundary  conditions  detined  in 


the  ldata 


file. 


OPTIONS:  None 

NOTE:  Descriptor  is  printed  as  an  identifier  at  the  beginning  of  the  odata  tile. 


UNE  2.  *LOAD  CONTROL,  #  steps  total,  #steps  to  execute,  analysis  type 
Number  of  fields:  3 
FORMAT:  (IX,  715) 

FIELD  1 :  Number  of  load  steps  in  the  analysis. 

OPTIONS:  If  FIELD  value  is  negative,  RAMP  function  is  flagged  ON. 

NOTE:  If  value  is  positive,  loads  (nodal,  pressure,  thermal,  etc.)  must  be  defined  below 
for  each  load  step.  If  value  is  negative,  only  the  total  of  each  load  type  is  read  and 
uniform  load  steps  are  automatically  generated. 

FIELD  2:  Number  of  load  steps  to  be  executed. 

OPTIONS:  None  ^  , 

NOTE:  Situations  may  arise  where  a  number  of  load  steps  are  defined  but  the  user  only 

wants  an  initial  portion  of  those  executed  to  check  the  progress  of  the  analysis. 

FIELD  3:  Type  of  analysis. 

OPTIONS: 

a)  0  =  Linear  elastic  analysis, 

b)  1  =  Reserved  for  future  use, 

c)  2  =  Data  check.  ' 

NOTE:  Option  c)  echos  all  user’s  input  files  to  odata  file  in  an  easy  to  read  and  verity 

form. 


UNE  3.  *LOAD  INDEX,  #BCs,  #N-force,  #N-displ,  #MPCs,  #E-Press,  B-force,  T-loads 
Number  of  fields:  7 

FORMAT:  (IX,  715)  j  .  . ,  .  ' 

FIELD  1:  Number  of  lines  of  homogeneous  nodal  displacement  data  to  be  read  below 

LINE  4. 

OPTIONS:  Value  can  be  zero. 

FIELD  2:  Number  of  lines  of  nodal  force  data  to  be  read  below  LINE  5. 

OPTIONS:  Value  can  be  zero. 

FIELD  3:  Number  of  lines  of  nonhomogeneous  nodal  displacement  data  to  be  read 
below  LINE  6. 

OPTIONS:  Value  can  be  zero.  TT1_, 

FIELD  4:  Number  of  lines  of  multi-point  constraint  data  to  be  read  below  LlNt  o. 

OPTIONS:  Value  can  be  zero. 

FIELD  5:  Number  of  lines  of  element  pressure  data  to  be  read  below  LINE  7. 

OPTIONS:  Value  can  be  zero. 

FIELD  6:  Number  of  lines  of  body  force  data  to  be  read  below  LINE  8. 

OPTIONS:  Value  can  be  zero. 

FIELD  7:  Number  of  lines  of  thermal  load  data  to  be  read  below  LINE  9. 

OPTIONS: 
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a)  0  =  No  thermal  load, 

b)  1  =  Uniform  temperature, 

c)  2  =  Prescribed  nodal  temperatures, 

LINE  4.  *FIXED  BC,  define  homogenous  nodal  dof  displacements. 

Number  of  fields:  2-5  depending  on  element  type. 

NOT^^Number  of  line’s  defining  nodes  with  prescribed  homogeneous  displacements 
must  agree  with  FIELD  1  in  LINE  3. 

FIELD  1:  Node  number. 

OPTIONS:  None  . . 

FIELD  2:  Constraint  code  for  translation  in  global  x  direction. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  displacement). 

FIELD  3:  Constraint  code  for  translation  in  global  y  direction. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  displacement). 

FIELD  4:  Constraint  code  for  translation  in  global  z  direction. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  displacement).  ,  ATT.nc, 

NOTE:  Not  required  for  two  dimensional  elements  (QUAD4  and  QUAD  ) 

FIELD  5:  Constraint  code  for  rotation  about  global  x-axis. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  rotation).  .  J  . 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 

FIELD  6:  Constraint  code  for  rotation  about  global  y  axis. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  rotation).  .  J  . 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 


LINE  5.  *NODAL  FORCES,  define  nodal  dof  forces. 

Number  of  fields:  2-5  depending  on  element  type. 

FORMAT:  (IX,  15,  5X,  5(G10.4))  .  .  n  ~ 

NOTE  1:  Number  of  lines  defining  nodes  with  prescribed  loads  must  agree  with  FIELD 

NOTE2Number  of  lines  defining  nodal  loads,  specified  by  FIELD  2  in  LINE  3  must  be 
listed  for  each  load  step  specified  in  FIELD  1  in  LINE  2  unless  RAMP  option  is  ON. 
FIELD  1:  Node  number. 

OPTIONS:  None 


FIELD  2:  Force  in  global  x  direction. 
OPTIONS:  Can  have  zero  value. 
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FIELD  3:  Force  in  global  y  direction. 

OPTIONS:  Can  have  zero  value. 

FIELD  4:  Force  in  global  z  direction. 

OPTIONS:  Can  have  zero  value.  ,^TTATAyl  .nilATV1, 

NOTE:  Not  required  for  two  dimensional  elements  (QUAD4  and  QUADS). 
FIELD  5:  Moment  about  global  x-axis. 

OPTIONS:  Can  have  zero  value.  .  j  , 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 
FIELD  6:  Moment  about  global  y-axis. 

OPTIONS:  Can  have  zero  value.  , 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 


LINE  6.  *NODAL  DISPS,  define  nonzero  nodal  displacements. 

Number  of  fields:  4-1 1  depending  on  element  type. 

FORMAT:  (IX,  15,  IX,  A5,  5(G10.4)) 

NOTE  1 :  Number  of  lines  defining  nodes  with  prescribed  nonzero  displacements  must 
agree  with  FIELD  3  in  LINE  3  unless  RAMP  function  is  flagged  ON. 

NOTE  2:  The  same  number  of  lines  defining  nodal  displacements,  specified  by 
in  LINE  3,  must  be  listed  for  each  load  step  specified  in  FIELD  1  in  LINE  2  unless 

RAMP  option  is  ON.  _  . ,  ., 

NOTE  3:  Prescribed  displacement,  both  zero  and  nonzero,  defined  here  ovemde  any 

defined  for  the  same  node  number  in  LINE  4. 

NOTE  4:  FIELDS  2-6  are  read  as  a  single  word  and  converted  to  integers  internally. 
FIELD  1:  Node  number. 


OPTIONS:  None 

FIELD  2:  Constraint  code  for  translation  in  global  x  direction. 


OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  displacement), 

c)  2  =  nonzero  displacement. 

FIELD  3:  Constraint  code  for  translation  in  global  y  direction. 
OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  displacement), 

c)  2  =  nonzero  displacement. 

FIELD  4:  Constraint  code  for  translation  in  global  z  direction. 


OPTIONS: 

a)  0  =  unconstrained 

b)  1  =  fixed  (zero  displacement) 

c)  2  =  nonzero  displacement  T  . 

NOTE:  Not  applicable  for  two  dimensional  elements  (QUAD4  and  QUAD8) 
FIELD  5:  Constraint  code  for  rotation  about  global  x-axis. 


OPTIONS: 

a)  0  =  unconstrained 

b)  1  =  fixed  (zero  rotation) 

c)  2  =  nonzero  rotation 
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NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 

FIELD  6:  Constraint  code  for  rotation  about  global  y-axis. 

OPTIONS: 

a)  0  =  unconstrained, 

b)  1  =  fixed  (zero  rotation), 

c)  2  =  nonzero  rotation.  . 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 

FIELD  7:  Translation  in  global  x  direction. 

OPTIONS:  Can  have  zero  value. 

FIELD  8:  Translation  in  global  y  direction. 

OPTIONS :  Can  have  zero  value. 

FIELD  9:  Translation  in  global  z  direction. 

OPTIONS:  Can  have  zero  value.  _ 

NOTE:  Not  applicable  for  two  dimensional  elements  (QUAD4  and  QUAD8) 

FIELD  10:  Rotation  about  global  x  axis. 

OPTIONS:  Can  have  zero  value.  . 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 

FIELD  11:  Rotation  about  global  y  axis. 

OPTIONS:  Can  have  zero  value. 

NOTE:  Currently  only  used  with  SHELL8  elements,  not  required  otherwise. 

LINE  7.  *MPC,  define  multi-point  dof  constraints. 

Number  of  fields:  8-29  depending  on  the  constraint  equation  formulated. 

FORMAT:  (IX,  15,  G12.5,3(I5,I2,G12.5),:,2(/,1 8X, 3(15,12, G12.5)) 

NOTE  1:  Number  of  lines  defining  multi-nodal  degree-of-freedom  constraint 
relationships  must  agree  with  FIELD  4  in  LINE  3. 

NOTE  2:  Up  to  three  lines  may  be  required  to  define  a  single  constraint  equation. 

NOTE  3-  Constraint  equations  allow  the  user  to  define  complex  relationships  between 
degrees-of-freedom  of  different  nodes.  For  example,  constrain  equations  are  useful  in 
applying  boundary  conditions  representing  a  plane  of  symmetry  of  a  model.  Constraint 
equations  have  the  form: 

C,  *DOF?  +C2  * DOFj  +...  +  C„  *DOFkn  =  C 
where : 

C,  =  Constants 

DOF"  =  Displacement  in  degree  of  freedom  i  of  node  n 

A  simple  constraint  equation  coupling  the  degrees-of-freedom,  in  the  x  direction  of  nodes 
1  and  2,  so  they  move  together  (same  magnitude),  would  look  like: 

(l)  *  DOF'  +  (- 1) *  DOF?  =  0 

FIELD  1:  Number  of  terms  in  constraint  equation. 

OPTIONS:  2  (min)  -  9  (max) 

FIELD  2:  Right  hand  value  (constant)  of  constraint  equation. 
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OPTIONS:  None. 

FIELD  3:  Node  number  of  first  term  in  constraint  equation. 
OPTIONS:  None 

FIELD  4:  Degree-of- freedom  constrained. 

OPTIONS:  1-5  depending  on  element  type. 

FIELD  5:  Coefficient  of  first  term  in  constraint  equation. 


OPTIONS:  None 

FIELDS  6-ll:Repeat  FIELDS  3-5  for  the  next  two  terms  in  the  constraint  equation. 
OPTIONS:  None 

NOTE  1 :  FIELDS  6-11  are  on  the  same  line  as  FIELDS  1-5. 

NOTE  2:  Number  of  terms  must  agree  with  FIELD  1 

FIELDS  12-20:  Repeat  FIELDS  3-5  for  the  next  three  terms  in  the  constraint  equation. 
OPTIONS:  None 

NOTE  1:  FIELDS  12-20  are  all  on  a  second  line  below  FIELDS  1-11. 

NOTE  2:  Number  of  terms  must  agree  with  FIELD  1 

FIELDS  21-29:  Repeat  FIELDS  3-5  for  the  next  three  terms  in  the  constraint  equation. 
OPTIONS:  None 

NOTE  1:  FIELDS  21-29  are  all  on  a  third  line  below  FIELDS  12-20. 

NOTE  2:  Number  of  terms  must  agree  with  FIELD  1 


LINE  8.  *ELM  PRESSURES,  define  element  edge  or  face  pressures. 

Number  of  fields:  4-7  depending  on  element  type. 

FORMAT:  (IX,  215, 4G12.4)  ‘  roTn, 

NOTE  1:  Number  of  lines  defining  applied  element  pressures  must  agree  with  FIELD  5 

NOTE  2:  Number  of  lines  defining  element  pressures,  specified  by  FIELD  5  in  LINE  3, 
must  be  listed  for  each  load  step  specified  in  FIELD  1  in  LINE  2  unless  RAMP  option  is 

ON. 

NOTE  3:  Element  pressure  are  defined  as  positive  acting  into  an  element  face  or  edge. 
NOTE  4:  Element  faces  or  edges  are  defined  by  node  number  sequence  determined  by 
element  connectivity.  Apply  right-hand  rule  to  node  sequences  listed  in  Table  F.2 
through  Table  F.4  to  determine  the  NEGATIVE  normal  (positive  pressure)  to  edge  or 

NOTE  5:  Pressures  can  only  be  specified  at  comer  nodes.  Pressures  at  midside  nodes  (if 
any)  are  interpolated  internally.  Pressures  can  vary  linearly  over  an  element  face  or  edge. 
NOTE  6:  For  plane  elements  (QUAD4,  QUAD8,)  node  numbering  is  shown  in  Fig.  F.5. 
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The  element  edges  for  plane  elements  are  defined  in  Table  F.2 

Table  F.2  Edge  numbering  for  QUAD4  and  QUAD8  elements. 


Edge’ 

Comer  Nodes 

1 

U 

2 

J,K 

3 

r  K.L 

4 

L,I 

Element  edges  and  faces  for  SHELL8  element  (Fig.  F.5)  are  defined  in  Table  F.3. 


Table  F.3  Edge  and  face  numbering  for  SHELL8  element. 


Edge  or  Face 

Comer  Nodes 

1 

I,J,K,L 

2 

U 

3 

J,K 

4 

K,L 

5 

L,I 

6 

I,L,K,J 

For  BRICK8  element  node  sequence  is  shown  in  Fig.  F.6. 
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Fig.  F.6  Node  numbering  for  BRICK8  element. 


The  element  faces  for  BRICK8  element  are  defined  in  Table  F.4. 


Table  F.4  Face  numbering  for  BRICK8  element. 


Edge  or  Face 

Comer  Nodes 

1 

I,J,K,L 

2 

I,M,N,J 

3 

J,N,0,K 

4 

K,0,P,L 

5 

I.L,P,M 

6 

M,P,0,N 

FIELD  1:  Element  number. 

OPTIONS:  None 

FIELD  2:  Face  or  edge  number. 

OPTIONS:  None. 

NOTE:  See  NOTE  6  for  face  and  edge  numbering. 
FIELD  3:  Element  number. 

OPTIONS:  None 

FIELD  4:  Pressure  at  first  node  defining  face  or  edge. 
OPTIONS:  None. 

FIELD  5:  Pressure  at  second  node  defining  face  or  edge. 
OPTIONS:  None. 

Note:  Only  required  for  specifying  face  pressures. 
FIELD  6:  Pressure  at  third  node  defining  face. 
OPTIONS:  None. 

Note:  Only  required  for  specifying  face  pressures. 
FIELD  7:  Pressure  at  fourth  node  defining  face  or  edge. 
OPTIONS:  None. 

NOTE:  Only  required  for  specifying  face  pressures. 
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LINE  9.  *BODY  FORCE,  define  forces  due  to  accelerations. 

Number  of  fields:  5-11  depending  on  acceleration  type. 

FORMAT:  (IX,  3G10.4)  .LI_Tn,.  T  tmtt 

NOTE  1:  The  existence  of  lines  defining  body  forces  must  agree  with  FIELD  6  in  LIN 


NOTE  2'  If  both  rotational  and  translational  body  forces  are  defined,  FIELDS  3-8  are 
used  to  define  rotation  vectors  and  FIELDS  9-1 1  are  used  to  define  direction  vector. 
NOTE  3:  Body  forces  can  only  be  applied  in  an  analyses  using  a  single  load  step. 
NOTE  4:  A  density  must  be  defined  for  each  material  in  mdata. 

FIELD  1:  Rotational  acceleration. 

OPTIONS :  Can  have  zero  value. 

FIELD  2:  Translational  acceleration. 

OPTIONS:  Can  have  zero  value. 

NOTE  1:  This  value  is  generally  gravity  or  some  multiple  of  it. 

NOTE  2:  FIELDS  1-2  form  the  first  of  2-4  lines  of  input  defining  body  forces. 

FIELD  3:  Global  x-coordinate  of  a  vector. 

OPTIONS: 

a)  Free  vector  tail  if  defining  axis  of  rotation, 

b)  Bound  vector  head  if  defining  direction  of  translation. 

NOTE:  OPTION  b)  is  used  if  only  translational  acceleration  exists. 

FIELD  4:  Global  y-coordinate  of  vector. 

OPTIONS 

a)  Free  vector  tail  if  defining  axis  of  rotation, 

b)  Bound  vector  head  if  defining  direction  of  translation. 

NOTE:  OPTION  b)  is  used  if  only  translational  acceleration  exists. 

FIELD  5:  Global  z-coordinate  of  vector. 

OPTIONS: 

a)  Free  vector  tail  if  defining  axis  of  rotation, 

b)  Bound  vector  head  if  defining  direction  of  translation. 

NOTE:  OPTION  b)  is  used  if  only  translational  acceleration  exists. 

NOTE:  FIELDS  3-5  form  the  second  group  of  2-4  lines  of  input  defining  body  forces. 
FIELD  6:  Global  x-coordinate  of  free  vector  head  that  defining  axis  of  rotation. 
OPTIONS:  None. 

NOTE:  Required  for  rotational  body  force  only. 

FIELD  7:  Global  y-coordinate  of  free  vector  head  defining  axis  of  rotation. 
OPTIONS:  None. 

NOTE:  Required  for  rotational  body  force  only. 

FIELD  8:  Global  z-coordinate  of  free  vector  head  defining  axis  of  rotation. 
OPTIONS:  None. 

NOTE:  Required  for  rotational  body  force  only. 

NOTE:  FIELDS  6-8  form  the  third  group  of  2-4  lines  of  input  defining  body  forces. 
FIELD  9:  Global  x-coordinate  of  bound  vector  defining  direction  of  translation. 
OPTIONS:  None. 

NOTE:  Required  only  if  both  rotational  and  translational  body  forces  exist. 

FIELD  10:  global  y-coordinate  of  bound  vector  defining  direction  of  translation. 
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OPTIONS:  None.  .  .  r 

NOTE:  Required  only  if  both  rotational  and  translational  body  forces  exist. 

FIELD  11:  Global  z-coordinate  of  bound  vector  defining  direction  of  translation. 

OPTIONS:  None.  ,  L  .  . 

NOTE  1 :  Required  only  if  both  rotational  and  translational  body  forces  exist. 

NOTE  2:  FIELDS  9-1 1  form  the  fourth  group  of  2-4  lines  of  input  defining  body  forces. 


LINE  10.  *TEMPERATURE  LOAD,  define  uniform  or  nodal  temperatures 
Number  of  fields:  2-3  depending  on  option  requested. 

FORMAT:  Either  (IX,  5G10.4)  (uniform  temperature  option)  or  (IX,  15, 2G1U.5) 


(prescribed  nodal  temperatures  option) 

NOTE  1 :  Number  of  lines  defining  thermal  loads  must  agree 


with  FIELD  6  in  LINE 


3. 


NOTE  2'  If  uniform  temperature  option  is  requested,  only  one  line  of  input  required.  If 
nodal  temperature  option  is  requested,  temperature  at  every  node  must  be  specified. 
NOTE  3:  Number  of  lines  defining  thermal  loads,  specified  by  FIELD  6  m  LINE  3,  must 
be  listed  for  each  load  step  specified  in  FIELD  1  in  LINE  2  unless  RAMP  option  is  ON. 
NOTE  4:  Analysis  requiring  iterations  for  load  equilibrium  (failure  or  nonlinear  shear), 
all  thermal  loads  are  applied  in  the  first  load  step. 

FIELD  1 :  Either  temperature  (uniform  option)  or  node  number  (nodal  option). 


OPTIONS:  None.  .  J  ,  .  .  ,  , 

FIELD  2:  Either  reference  temperature  (uniform  option)  or  nodal  temperature  (no  a 


option). 

OPTIONS:  None. 

FIELD  3 :  Reference  temperature  (nodal  option). 
OPTIONS:  None. 

NOTE:  Not  required  under  the  uniform  temperature  option. 


LINE  11.  *OUTPUT  CONTROL,  Maxpm,  pmsteps,  Local,  Gauss,  Lay,  Fail,  ANSYS, 
Plot 

Number  of  fields:  8. 

FORMAT:  (IX,  915)  .  ^  ,  .  u  .  ,  , . 

NOTE  1:  First  of  two  LINES  controlling  form  and  quantity  of  analysis  results  printed  in 

od-sitsi 

NOTE  2:  Generally,  descriptors  of  stress  and  strain  values,  based  on  flags  set  in  LINE 

1 1 ,  are  printed  in  odata. 

FIELD  1:  Request  gross  printing  of  results. 

OPTIONS: 

a)  0  =  Gross  print  flag  is  OFF,  m  n 

b)  1  =  All  load  steps  are  printed  using  print  control  pattern  requested  in  LINE  12, 

c)  >1  =  Load  steps  that  are  multiples  of  this  value  will  be  printed  using  print  control 
pattern  requested  in  LINE  12. 

FIELD  2:  Number  of  load  steps  specifically  flagged  to  be  printed  using  print  contro 
patterns  requested  for  each  load  step  listed  in  LINE  12. 

OPTIONS: 

a)  0  =  Load  step  print  flag  is  OFF, 
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b)  >1  =  Number  of  load  steps  to  be  printed.  Load  step  print  flag  is  ON. 
NOTE:  FIELD  2  is  over-ridden  by  FIELD  1  if  gross  print  flag  is  ON. 
FIELD  3:  Coordinate  system  reference  for  stress  and  strain  output. 


OPTIONS: 

a)  0  =  Global  coordinate  system, 

b)  l  =  Principal  material  (fiber  direction)  coordinate  system. 
FIELD  4:  Location  within  an  element  for  stress  and  strain  output. 

OPTIONS: 

a)  0  —  Nodes, 

b)  1  =  Gauss  points. 

FIELD  5:  Location  within  a  layer  for  stress  and  strain  output. 


OPTIONS: 

a)  -1  =  Bottom  of  a  layer, 

b)  0  =  Midplane, 

c)  1  =  Top  of  a  layer. 

FIELD  6:  Request  failure  analysis. 


OPTIONS: 

a)  0  =  Failure  analysis  flagged  OFF. 

b)  l  =  Predict  an  ultimate  failure  load  for  each  element.  See  NOTES  1  &  . 

c)  2  =  Fail  fiber  and  matrix  according  to  criteria  set  in  fdata.  See  NOTE  3 . 

d)  3  =  Predict  composite  failure  according  to  3D  Tsai-Wu  criterion.  See  NOTE  3. 
NOTE  1 :  OPTION  b)  requires: 

a)  Single  load  step, 

b)  One  load  type  per  analysis,  e.g.,  nodal  forces, 

c)  Magnitude  of  applied  load  must  equal  1  (unit  load), 

d)  QUAD4  or  QUAD8  elements  only, 

e)  ANSYS  output  flag  must  be  set  ON  (see  FIELD  7). 

NOTE  2:  Ultimate  failure  loads  are  output  as  ETABLE  entries  in  the  ANSYS  readable 

output  files.  See  FIELD  7.  '  ,  ,  «  ■ 

NOTE  3:  OPTION  b)  requires  fdata  file  and  incremental  (multiple  load  steps)  analysis. 

Can  be  used  with  any  element  type. 


FIELD  7:  Request  output  files  of  stress,  strain  and  displacement  results  in  an  ANSYS 

readable  format. 

OPTIONS: 


a)  0  =  ANSYS  output  flagged  OFF. 

b)  >0  =  ANSYS  output  flagged  ON.  See  NOTES  1  and  2. 

NOTE  1 :  OPTION  b)  causes  three  or  four  output  files  to  be  created: 

1  .mct2ans.c  (composite  displacements,  strains  and  stresses), 

2.mct2ans.r  (reinforcement  displacements,  strains  and  stresses), 

3  mct2ans.m  (matrix  displacements,  strains  and  stresses). 

4  damage.dat  (material  number  corresponding  to  composite  damage  state).  _ 

If  ultimate  failure  loads  are  requested  (FIELD  6  OPTION  C)  they  are  output  as  ETABLE 

entries  in  each  of  the  three  files.  ,  ...  .  ., 

NOTE  2-  In  a  multi-step  analysis  MCT_77  writes  ANSYS  readable  out  put  starting  a 
load  step  specified,  e.g.,  a  value  of  3  will  output  all  results  from  load  step  3  onward. 
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FIELD  8:  Request  output  file  of  results  specified  in  LINE  13  in  a  format  suitable  to  be 
read  by  graphing  program  or  spreadsheet. 

OPTIONS: 

a)  0  =  Plot  output  flagged  OFF. 

b)  1  =  Plot  output  flagged  ON.  If  stress  or  strain  are  requested,  composite  values  will 

c)  2=Upiot  output  flagged  ON.  If  stress  or  strain  are  requested,  reinforcement  values 

will  be  used.  ,  .n  u 

d)  3  =  Plot  output  flagged  ON.  If  stress  or  strain  are  requested,  matrix  values  will  be 

used. 

NOTE  :  OPTIONS  b)-d)  create  an  ASCII  format  file,  plot.dat. 


LINE  12.  *PRINT  CONTROL,  Inc,  Displs,  Nsts,  Nstn,  Asts,  Astn,  Ests,  Estn,  Psts,  Ists, 
Rfor 

Number  of  fields:  1 1 . 

FORMAT:  (IX,  15,  IX,  A10)  .  r  ,  .  u  ... 

NOTE  1:  Second  of  two  LINES  controlling  form  and  quantity  of  analysis  results  printed 

NOTE  2:  Input  under  LINE  13  is  required  only  if  either  FIELD  1  or  2  m  LINE  1 1  is 

NOTlf  3^ Number  of  lines  defining  print  control  patterns  must  agree  with  FIELDS  1-2  in 


NOTE4:‘  In  most  cases,  setting  flags  specified  in  LINE  12  to  OFF  causes  program  to  skip 
execution  of  associated  subroutines,  significantly  decreasing  overall  execution  times  for 


large  models. 

FIELD  1 :  Load  step  to  be  printed 
OPTIONS:  None. 

NOTE  1:  Ignored  if  gross  print  flag  is  set  ON  in  FIELD  1  of  LINE  11. 
FIELD  2:  Print  request  -  displacements. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes). 

FIELD  3:  Print  request  -  nodal  stresses. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Y es). 

FIELD  4:  Print  request  -  nodal  strains. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes). 

FIELD  5:  Print  request  -  stresses  averaged  at  nodes. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes). 

FIELD  6:  Print  request  -  strains  averaged  at  nodes. 

OPTIONS: 
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a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Y es). 

FIELD  7:  Print  request  -  stresses  calculated  at  element  centroid. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Y es). 

FIELD  8:  Print  request  -  strains  calculated  at  element  centroid. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Y es). 

FIELD  9:  Print  request  -  principal  stresses. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes).  ... 

FIELD  10:  Print  request  -  transversely  isotropic  stress  invariants. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes). 

FIELD  11 :  Print  request  -  reaction  forces  at  nodes. 

OPTIONS: 

a)  N  =  Print  request  is  OFF  (No), 

b)  Y  =  Print  request  is  ON  (Yes). 


LINE  13.  *PLOT  CONTROL,  {load  id,  load#,  blank,  blank}  or  {field  id,  node,  lay, 
elem} 

Number  of  fields:  8. 

FORMAT:  (IX,  2(415))  .  ITl^11.'  '  JAM 

NOTE  1 :  Input  under  LINE  13  is  required  only  if  FIELD  8  in  LINE  1 1  is  flagged  ON. 
NOTE  2:  Typically  two  corresponding  quantities  are  output,  e.g.,  stress  and  strain. 
NOTE  3:  LINE  13  specifies  incremental  quantities  to  be  written.  Quantities  that  can  e 


written  are:  . ,  ,  ,  , 

a)  Field  variables  LINE  12  using  FIELDS:  {field  identifier,  node  number,  layer, 

b)  LoSshi  LINE  3  using  FIELDS:  {load  identifier,  load  number,  BLANK,  BLANK} 

The  field  identifiers  are: 

a)  1  =  Nodal  displacements, 

b)  2  =  Nodal  stresses, 

c)  3  =  Nodal  strains, 

d)  4  =  Stresses  averaged  at  nodes, 

e)  5  =  Strains  averaged  at  nodes, 

f) 6  =  Stresses  at  element  centroid, 

g)  7  =  Strains  at  element  centroid, 

h)  8  =  Principal  stresses, 

i)  9  =  Transversely  isotropic  stress  invariants. 

All  field  variables  of  that  type  and  at  that  point  are  printed,  e.g.  displacements  in  the  x,  y, 

and  z  direction. 
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The  load  identifiers  are: 

a)  21  =  Nodal  loads, 

b)  22  =  Prescribed  (nonzero)  nodal  displacements, 

c)  23  =  Element  pressures, 

d)  24  =  Nodal  temperature. 

Load  number  corresponds  to  position,  within  loads  of  that  type,  as  specified  in  ldata 
e.g.,  first  element  pressure,  second  element  pressure,  tenth  element  pressure,  etc. 

FIELD  1 :  Quantity  identifier. 

OPTIONS: 

a)  Field  identifier, 

b)  Load  identifier. 

NOTE:  See  NOTE  3  under  LINE  13. 

FIELD  2:  Sub-identification. 

OPTIONS: 

a)  Node  number  (field  option), 

b)  Load  number  (load  option). 

FIELD  3:  Sub-identification. 

OPTIONS: 

a)  Layer  number  (field  option), 

b)  Blank  (load  option). 

FIELD  4:  Sub-identification. 

OPTIONS: 

a)  Element  number  (field  option), 

b)  Blank  (load  option). 

FIELDS  5-8:  Repeat  FIELDS  1-4  for  a  corresponding  quantity. 

OPTIONS:  Same  as  for  FIELDS  1-4. 

NOTE:  Fields  and  loads  can  be  listed  as  corresponding  quantities. 
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F.6.  Mdata  File 

The  mdata  file  contains  information  defining  elastic  material  constants  for  up  to  eight 
orthotropic  materials.  A  sample  mdata  file  is  shown  in  Fig.  F.7.  The  mdata  file  requires 
a  minimum  of  three  LINES  to  define  a  single  homogeneous  material,  a  minimum  of  5 
LINES  to  define  a  composite  material  and  its  two  constituents,  and  a  minimum  of  9  lines 
to  define  a  single  composite  material  and  its  two  constituents  in  damaged  and  undamaged 
states.  All  materials  are  assumed  to  be  transversely  isotropic.  A  LINE-by-LINE 
description  follows. 


LINE  1. '♦Material  Title 
Number  of  fields:  1 

FORMAT:  (A60)  . 

FIELD  1 :  General  description  of  the  materials  defined  in  the  mdata  file. 

OPTIONS:  None  .  t 

NOTE:  Descriptor  is  printed  as  an  identifier  at  the  beginning  of  the  odata  tile. 

LINE  2.  *Material  Sets 
Number  of  fields:  2 
FORMAT:  (IX,  915) 

FIELD  1:  Number  of  material  constant  sets  to  be  read. 

OPTIONS:  1  (min)  -  8  (max). 

FIELD  2:  Request  analysis  with  nonlinear  shear  effects. 

OPTIONS: 

a)  0  =  Nonlinear  shear  analysis  flagged  OFF, 

b)  1  =  Nonlinear  shear  analysis  flagged  ON.  See  NOTE. 

NOTE:  OPTION  b)  requires  an  sdata  file  and  incremental  (multiple  load  step)  analysis. 
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IVIdici  idi  i  iiiw  . 

Boron/5505,  50%  FVF,  Order  important  for  failure  analysis 
•Material  Sets 
7  1 

•Material, 1  .Undamaged  composite 
4  5  0.50 

207.0D+9  20.98D+9  20.98D+9 
0  29369  0.29369  0.45140 

7  715D+9  7.71 5D+9  7.229D+9 
0.1012D-4  0.2803D-4  0.2803D-4 
2.000D+3  O.OOOD+O  O.000D+0 
•Material, 2,Composite  with  damaged  matrix 
4  7  0.50 

203.3D+9  0.2201D+9  0.2201  D+9 
0.29349  0.29349  0.46414 

0  8050D+8  0.8050D+8  0.751 5D+8 
0.9526D-5  0.2812D-4  0.2812D-4 
2.000D+3  O.OOOD-O  O.OOOD-O 
•Material, 3 .Composite  with  damaged  fiber 
6  5  0.50 

3  644D+9  1.733D+9  1.733D+9 

0.34846  0.34846  0.38027 

0  9089D+9  0.9089D+9  0.6282D+9 
0  4109D-4  0.4063D-4  0.4063D-4 
2.000D+3  O.OOOD-O  O.OOOD-O 
•Material, 4,Undamaged  fiber 

4  4  0.00 

407.0D+9  407.0D+9  407.0D+9 
0.25000  0.25000  0.25000 

162.8D+9  162.8D+9  162.8D+9 

0.9520D-5  0.9520D-5  0.9520D-5 
2.600D+3  O.OOOD-O  O.OOOD-O 
•Material.5 .Undamaged  matrix 

5  5  0.00 

7.239D+9  7.239D+9  7.239D+9 
0.35000  0.35000  0.35000 

2.681  D+9  2.681D+9  2.681D+9 
0.4120D-4  0.4120D-4  0.4120D-4 
1.200D+3  O.OOOD-O  O.OOOD-O 
•Material, 6, Damaged  fiber 

6  6  0.00 

407.0D+5  407.0D+5  407.0D+5 
0.25000  0.25000  0.25000 

162.8D+5  162.8D+5  162.8D+5 

0.9520D-5  0.9520D-5  0.9520D-5 
2.600D+3  O.OOOD-O  O.OOOD-O 
•Material, 7,Damaged  matrix 
7  7  0.00 

7  239D+7  7.239D+7  7.239D+7 

0.35000  0.35000  0.35000 

2.681  D+7  2.681  D+7  2.681  D+7 

0.4120D-4  0.4120D-4  0.4120D-4 
1.200D+3  O.OOOD-O  O.OOOD-O _ _ _ 


Fig.  F.7  Sample  mdatafile. 


LINE  3.  ^Material  # 

FOmT:  (IX,  215, 14G5.1)  for  FIELDS  1-3  and  (1X,6G12.7)  for  FIELDS  4-18. 
NOTE  1:  Six  lines  of  input  follow  LINE  3. 

NOTE  2:  Material  number  is  determined  by  its  sequential  position  within  Ida  a  i  e, 
NOTE  3:  The  existence  of  LINES  3-10  must  agree  with  FIELD  1  in  LINE  2. 
FIELD  1:  Material  number  of  corresponding  reinforcement  constituent. 


NOTE  •  All  materials  are  assumed  to  be  composites.  FIELD  1  specifies  which  material, 
listed  sequentially  within  ldata  file,  defines  reinforcement  material  constants  for  this 
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composite.  Homogeneous  materials  are  defined  by  specifying  constituent  material 
numbers  (FIELDS  1  and  2)  to  be  the  same  as  the  “composite”. 

FIELD  2:  Material  number  of  corresponding  matrix  constituent. 

NOTE  :  All  materials  are  assumed  to  be  composites.  FIELD  2  specifies  which  materia , 
listed  sequentially  within  ldata  file,  defines  matrix  material  constants  f°r  thls  composite.. 
Homogeneous  materials  are  defined  by  setting  constituent  material  numbers  (FIELDS  1 
and  2)  to  be  the  same  as  the  “composite”. 

FIELD  3:  Reinforcement  volume  fraction. 

OPTION:  None.  . 

NOTE  1 :  FIELD  3  for  homogenous  material  defaults  internally  to  u.du. 

NOTE  2:  FIELDS  1-3  form  the  first  of  6  lines  of  input  under  LINE  3. 

FIELD  4:  Young’s  modulus  along  longitudinal  axis  of  fiber  (En). 

OPTION:  None.  .  ,  . 

FIELD  5:  Young’s  modulus  transverse  to  longitudinal  axis  of  fiber  (£22)- 

OPTION:  None.  . 

FIELD  6:  Young’s  modulus  transverse  to  longitudinal  axis  of  fiber  (E33L 

OPTION:  None. 

NOTE  1 :  Under  the  assumption  of  transverse  isotropy,  E33  -£22- 
NOTE  2:  FIELDS  4-6  form  the  second  of  6  lines  of  input  under  LINE  3. 

FIELD  7:  Major  Poisson’s  ratio  (v]2). 

OPTION:  None. 

NOTE:  Defined  as: 


£2 

V"='T 


FIELD  8:  Major  Poisson’s  ratio  (V13). 

OPTION:  None. 

NOTE:  Under  the  assumption  of  transverse  isotropy,  vJ3  =vi2- 
FIELD  9:  Transverse  Poisson’s  ratio  (v23>. 

OPTION:  None. 

NOTE:  Defined  as: 


NOTE:  FIELDS  7-9  form  the  third  of  6  lines  of  input  under  LINE  3. 
FIELD  10:  Longitudinal  shear  modulus  (G12). 

OPTION:  None. 

FIELD  11:  Longitudinal  shear  modulus  (G13). 

OPTION:  None. 

NOTE  1 :  Under  the  assumption  of  transverse  isotropy,  G13-G12. 

FIELD  12:  Transverse  shear  modulus  (G23>- 

OPTION:  None.  t\tc  1 

NOTE:  FIELDS  10-12  form  the  fourth  of  6  lines  of  input  under  LINE  3. 

FIELD  13:  Coefficient  of  thermal  expansion  along  fiber  axis  (ai). 

OPTION:  None. 

FIELD  14  Coefficient  of  thermal  expansion  transverse  to  fiber  axis  (0,2). 
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OPTION:  None.  . 

NOTE  1 :  Under  the  assumption  of  transverse  isotropy,  Gb  Gb- 
FIELD  15:  Coefficient  of  thermal  expansion  transverse  to  fiber  axis  (0.3). 
OPTION:  None. 

NOTE  1 :  Under  the  assumption  of  transverse  isotropy,  0.3-012. 

NOTE  2:  FIELDS  13-15  form  the  fifth  of  6  lines  of  input  under  LINE  3. 

FIELD  16:  Material  density  (p) 

OPTION:  None. 

FIELD  17  BLANK 
OPTION:  None. 

FIELD  18:  BLANK 

OPTION:  None.  TTvn:i 

NOTE:  FIELDS  16-18  form  the  sixth  of  6  lines  of  input  under  LINE  3. 

LINE  4.  *Material  # 

Number  of  fields:  18  •  T  nc  ,  18 

FORMAT:  (IX,  215, 14G5.1)  for  FIELDS  1-3  and  (1X.6G12.7)  for  FIELDS  4-18. 

NOTE  1 :  LINES  4-10  are  identical  to  LINE  3,  each  defining  another  material. 
NOTE  2:  The  existence  of  LINES  4-10  must  agree  with  FIELD  1  in  LINE  2. 
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F.  7.  Sdata  File 

The  sdata  file  contains  parameters  of  curve  fit  equations  that  form  constitutive 
relationships  among  elastic  material  constants  for  a  composite  and  its  constituents  under 
nonlinear  shear  effects.  In  MCT_77,  only  one  material  can  have  nonlinear  shear 
behavior.  The  parameters  in  this  file  correspond  to  material  #1  in  the  mdata  file. 

A  changing  composite  shear  modulus  requires  incremental  application  of  the  load  and 
iterations  within  each  load  step  to  achieve  equilibrium.  Two  tolerances,  deflection 
(DTOL)  and  load  (FTOL)  control  equilibrium  iterations.  Default  values  (DTOL=7.0  and 
FTOL=1.0E-5)  are  used  in  a  non-failure  MCT_77  analysis  with  nonlinear  shear.  The 
user  can  change  the  tolerances  in  a  non-failure  analysis  by  setting  the  failure  analysis  flag 
to  ON  in  ldata  and  the  number  of  materials  to  0  in  fdata.  MCT_77  will  open  fdata  and 
read  user  defined  tolerance  values. 

A  sample  sdata  file  is  shown  in  Fig.  F.8.  The  sdata  file  requires  eleven  LINES.  A  LINE- 
by-LBSTE  description  follows. 

LINE  1.  *Curve  Fit  Title 
Number  of  fields:  1 

FORMAT:  (A60)  J .  .  „ 

FIELD  1 :  General  description  of  material  defined  in  sdata  tile. 

OPTIONS:  None  ,  .  r 

NOTE:  Descriptor  is  printed  as  an  identifier  at  the  beginning  of  odata  h  . 
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‘Curve  hit  i  itle 
Boron/5505,  50%  FVF 
‘Max  strain,  shear_12,  shear_23 
0.0462D0  0.0462D0 

‘Shear  12  stress/strain  parameters:  a_0,a_1,a_2,c_0,c_1,c_J4 
0.000D+0  -5.757D+1  -5.755D+1 
1.340D+8  -1.327D+8  -1.253D+6 
‘E22/G12  parameters:  p_0,p_1,p_2 
4.270D-2  3.853D+0  -1.201  D+0 

‘V12/G12  parameters:  p_0,p_1,p_2 
4.861  D-2  -1.290D-2  2.374D-3 

‘G23/G12  parameters:  p_0,p_1,p_2 
5.497D-3  1.095D+0  -1.674D-1 

‘G/G12  parameters:  p_0,p_1,p_2 
1.234D-5  3.327D-1  1.474D-2 

‘Shear  23  stress/strain  parameters:  a_0,a_1  ,a_2,c_0,c_i  ,c_z 
0  000D+0  -4.523D+1  -8.818D+2 
1.386D+8  -1.369D+8  -1.726D+6 
‘E22/G23  parameters:  p_0,p_1,p_2 
4.270D-2  3.853D+0  -1.201  D+0 

‘V12/G23  parameters:  p_0,p_1,p_2 
4.861  D-2  -1.290D-2  2.374D-3 

‘G12/G23  parameters:  p_0,p_1,p_2 
5.497D-3  1.095D+0  -1.674D-1 

‘G/G23  parameters:  p_0,p_1,p_2 

1.234D-5  3.327D-1  1.474D-2 _ _ _ 


Fig.  F.8  Sample  sdata  file. 


LINE  2.  *Max  strain,  shear_12,  shear_23. 
Number  of  fields:  2 


FIELIMf  Maximum  shear  strain  to  be  used  in  nonlinear  composite  shear  12  stress-strain 
curve-fit  equation.  See  NOTE 

F^ELD^  Maximum  shear  strain  to  be  used  in  nonlinear  composite  shear  23  stress-strain 
curve- fit  equation.  See  NOTE. 

OPTIONS:  None.  .  „  , 

NOTE:  Shear  strain  in  a  general  lamina  can  exceed  the  domain  of  shear  stress-strain 

equations  resulting  in  poor  equilibrium  convergence.  Shear  moduli  calculated  at  the 

maximum  defined  shear  strain  is  used  when  lamina  shear  strains  exceed  maximum  va  u  . 


UNE  3.  *Shear  12  stress/strain  parameters:  a_0,a_l,a_2,c_0,c_l,c_2 
Number  of  fields:  6 
FORMAT:  (1X,6G12.7) 

NOTE:  Defines  nonlinear  shear  12  stress/strain  curve. 

FIELD  1:  BLANK. 

OPTIONS:  None 

FIELD  2:  Curve  fit  parameter  in  first  exponent  (ai). 
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OPTIONS:  None. 

FIELD  3:  Curve  fit  parameter  in  second  exponent  w). 


OPTIONS:  None.  t  .  IIWC, 

NOTE  2:  FIELDS  1-3  form  the  first  of  2  lines  of  input  under  LINE  3. 

FIELD  4:  Curve  fit  constant  forming  first  term  of  equation  (c0). 


OPTIONS:  None.  r 

FIELD  5:  Curve  fit  constant  coefficient  in  second  term  of  equation  (cij. 

OPTIONS:  None.  * 

FIELD  6:  Curve  fit  constant  coefficient  in  third  term  of  equation  (C2). 

OPTIONS:  None. 

NOTE  1 :  Three  term  shear  stress/strain  curve  fit  equation  has  the  lorm. 


r  =  C0  +  Cxe(a'r)  +  C2el 


(<>2r) 


where  C,  and  a,  are  curve  fit  parameters  and  y  is  engineering  shear -strain. 
NOTE  2:  FIELDS  4-6  form  the  second  of  2  lines  of  input  under  LINE  3. 


LINE  4.  *E22/G12  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 
FORMAT:  (1X,6G12.7) 

NOTE:  Defines  relationship  between  E22  and  nonlinear  G12. 

FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po). 

OPTIONS:  None.  r  '  .  ,  . 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (Pu- 

OPTIONS:  None.  .  .  . 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (P2j- 

OPTIONS:  None.  . 

NOTE  1 :  Second  order  polynomial  curve  fit  equation  has  the  form. 

Value  =  p0+px*  ratio  +  p2  *  ratio2 

LINE  5.  *V12/G12  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 
FORMAT:  (IX, 6G12.7) 

NOTE:  Defines  relationship  between  vJ2  and  nonlinear  Gi2- 
FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po). 

OPTIONS:  None.  ■  r  '  .  ,  . 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (Pu 

OPTIONS:  None.  .  .  .  . 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (P2). 

OPTIONS:  None. 

LINE  6.  *G23/G12  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 
FORMAT:  (IX, 6G12.7) 

NOTE:  Defines  relationship  between  G23  and  nonlinear  G12. 
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FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po)- 

OPTIONS:  None.  f  .  ,  , 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (Pi) 

OPTIONS:  None.  .  .  . 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (p2). 

OPTIONS:  None. 


LINE  7.  *G/G12  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 


FORMAT:  (1X.6G12.7)  ,  .. 

NOTE:  Defines  relationship  between  G  (matrix  shear  modulus)  and  nonli 

FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po). 

OPTIONS:  None.  .  .  ,  . 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (pi). 

OPTIONS:  None.  r  ,  .  ,  , 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  the  equation  (p2). 

OPTIONS:  None. 


Gi2- 


LINE  8.  *  Shear  23  stress/strain  parameters:  a_0,a_l,a_2,c_0,c_l,c_2 
Number  of  fields:  6 
FORMAT:  (1X.6G12.7) 

NOTE:  Defines  nonlinear  shear  23  stress/strain  curve. 

FIELD  1:  BLANK. 

OPTIONS:  None 

FIELD  2:  Curve  fit  parameter  in  first  exponent  (ai). 

OPTIONS:  None. 

FIELD  3:  Curve  fit  parameter  in  second  exponent  (a2). 

OPTIONS:  None.  .  T_„c 

NOTE  2:  FIELDS  1-3  form  the  first  of  2  lines  of  input  under  LINE  8. 
FIELD  4:  Curve  fit  constant  forming  first  term  of  equation  (co). 
OPTIONS:  None. 

FIELD  5:  Curve  fit  constant  coefficient  in  second  term  of  equation  (ci). 

OPTIONS:  None.  ...  *  .  .  . 

FIELD  6:  Curve  fit  constant  coefficient  in  third  term  of  equation  (C2). 

OPTIONS:  None.  .  .  .  . 

NOTE  1 :  The  three  term  shear  stress/strain  curve  fit  equation  has  the  lorm 


T  =  C0+C/°'r)  +  C2eM 


where  C,  and  a,  are  curve  fit  parameters  and  yis  engineering  shear  strain. 

NOTE  2:  FIELDS  4-6  form  the  second  of  2  lines  of  input  under  LINE  8. 

NOTE  3:  If  shear  23  data  is  unavailable,  shear  12  data  can  be  used  as  an  approximation. 

LINE  9.  *E22/G23  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 
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FORMAT:  (1X.6G12.7) 

NOTE:  Defines  relationship  between  E22  and  nonlinear  023- 
FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po)- 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (pi). 

OPTIONS:  None.  ,  .  ,  , 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (.P27- 

OPTIONS:  None.  ...... 

NOTE  :  Second  order  polynomial  curve  fit  equation  has  the  lorm. 

Value  =  p0+ Pi*  ratio  +  p2*  ratio 2 

LINE  10.  *V12/G23  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 
FORMAT:  (1X,6G12.7) 

NOTE:  Defines  relationship  between  V12  and  nonlinear  G23. 

FIELD  1 :  Curve  fit  constant  forming  first  term  of  equation  (po)- 

OPTIONS:  None.  .  , 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (Pit- 

OPTIONS:  None.  ^  . 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (P2F 

OPTIONS:  None. 

LINE  11.  *G12/G23  parameters:  p_0,p_l,p_2 
Number  of  fields:  3  j 
FORMAT:  (1X,6G12.7) 

NOTE:  Defines  relationship  between  G12  and  nonlinear  G23. 

FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po)- 

OPTIONS:  None.  .  . 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (Pit- 

OPTIONS:  None.  ^  . 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (P2). 

OPTIONS:  None. 

LINE  12.  *G/G23  parameters:  p_0,p_l,p_2 
Number  of  fields:  3 

FORMAT:  (IX, 6G12.7)  .  ,  r  MrC 

NOTE:  Defines  relationship  between  G  (matrix  shear  modulus)  and  nonlinear  023- 

FIELD  1:  Curve  fit  constant  forming  first  term  of  equation  (po)- 

OPTIONS: None.  .  > 

FIELD  2:  Curve  fit  constant  coefficient  in  second  term  of  equation  (pi). 

OPTIONS:  None.  .  .  , 

FIELD  3:  Curve  fit  constant  coefficient  in  third  term  of  equation  (P2). 

OPTIONS:  None. 
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appendix  g 

OVERVIEW  OF  MATERIAL  FAILURE  THEORIES 

G.l.  Categorizing  Failure  Analysis 

Failure  of  a  material  can  be  classified  by  yield  (ductile  material  failure),  or  fracture 
(brittle  material  failure).  These  failure  modes  can  occur  in  combinations  making  them 
sometimes  difficult  to  identify.  Furthermore,  a  material  can,  depending  on  temperature, 
stress  state,  or  rate  of  loading,  act  in  either  a  brittle  or  ductile  manner  [1].  But  over  a 
wide  range  of  conditions,  a  particular  material  typically  behaves  in  either  a  ductile  or 
brittle  manner.  Within  this  paper,  only  failure  under  quasi-static  load  conditions  is 

considered. 

G.2.  Brittle  Material  Failure  Criteria 
G.2.1.  Linear  Elastic  Fracture  Mechanics 

For  brittle  failure,  the  most  commonly  used  failure  criteria  are  maximum  stress  or  strain, 
and  linear  elastic  fracture  mechanics  (LEFM).  LEFM  attempts  to  determine  the  stress 
necessary  to  cause  a  pre-existing  through-thickness  crack  in  a  primarily  brittle  material  to 
grow  in  an  unstable  (self-propagating  or  catastrophic)  manner.  The  three  most  common 
approaches  to  LEFM  are  [2]  :  classical,  stress  intensity  factor,  and  strain  energy  release 
rate.  It  is  important  to  note  that  compressive  failure  of  a  material  is  not  addressed  by 

LEFM. 
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Classical  (Griffith  type)  LEFM  approach  to  a  mode  I  (opening)  crack  propagation  defines 
only  the  critical  applied  stress  and  has  the  form 

(G1) 

where 


cr  =  applied  stress, 
a  =  half-crack  length, 

E  =  Young’s  modulus, 
y  =  surface  energy  per  unit  area, 

with  the  critical  stress  defining  the  mode  I  fracture  toughness  as 


(G.2) 


The  stress  intensity  factor  LEFM  approach  develops  stress  values  about  the  crack  tip  and 


has  the  form: 


where 


K'  :/,W 
K‘  fM 
fM 


nr 

K, 


nr 


ft  =  trigonometric  functions  of  the  angle  6, 

Ki  =  (T-Jm  =  stress  intensity  factor  (mode  I). 


(G.3) 


In  application,  these  approaches  are  similar  in  that  the  classical  critical  fracture  toughness 
(or  the  critical  stress  intensity  factor  which  is  identical)  is  a  material  property  that  can  be 
determined  experimentally.  In  application,  once  the  fracture  toughness  is  known,  the 
critical  applied  stress  can  be  determined  based  on  an  experimentally  determined  critical 

crack  size. 
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Both  of  the  above  approaches  can  be  difficult  to  implement  for  complex  loading, 
geometry,  and  anisotropic  materials.  The  strain  energy  release  rate,  as  is  typical  of 
energy  approaches,  is  generally  easier  to  apply  to  a  broader  range  of  problems  and  proven 
to  be  a  powerful  tool  in  crack  growth  studies  [2].  The  strain  energy  release  rate  approach 
to  LEFM  has  the  form: 


where 


r 

'  It  da 


(G.4) 


P  =  load 

t  =  through  crack  thickness, 
s  =  ,  the  system  compliance, 

u  =  change  in  displacement  at  the  point  of  load  application  due  to 
crack  growth, 


and  the  critical  strain  energy  release  rate  for  mode  one  defined  as: 


2 1  \  da. 


The  advantage  of  the  strain  energy  release  rate  approach  is  that  knowledge  of  material 
properties  or  crack  tip  stresses  is  not  required.  Simply  plotting  system  compliance  verses 
crack  length  and  computing  the  slope  at  critical  load  determines  critical  strain  energy 


release  rate.  It  has  been  shown  [3]  that: 

K)c  =  GICE  (G.5) 

Hence,  the  critical  applied  stress  can  be  determined  based  on  a  known  critical  crack  size. 
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There  are  analogous  K,c  developments  for  craek  propagation  caused  by  mode  II  (sliding) 
and  mode  m  (tearing)  loads.  Load  and  boundary  conditions  can  produce  crack 
propagation  that  is  a  combination  of  cracking  modes  requiring  a  “mixed”  mode  analysis. 


G.2.2.  Maximum  Stress  (or  Strain)  Failure  Criteria 

The  maximum  stress  (strain)  criteria  can  be  stated  as  follows. 

“The  material  at  a  point  in  a  structure  fails  when  the  maximum  stress  (strain)  exceeds  the 
,e?sne  stess  (s.ra,n)  m  the  minimum  stress  (sham)  exceeds  the  maxmumt  compressrve 
stress  (strain)  as  determined  in  a  simple  uni-axial  specimen  oft  e  same  ma  era 
failed”  [4]. 


The  maximum  and  minimum  stresses  (strains)  stated  above  can  refer  to  principal 
(coordinate  invariant)  stresses  but  can  also  refer  to  local  stresses  in  principal  material 
directions  (coordinate  dependent)  typical  of  fiber  reinforced  composites. 

In  general,  a  material  will  not  fail  if. 

~S  <  <7,<+S 

~S  <  <t2<+5  < G  6) 

~S  <  a3<+S 


+S  =  ultimate  unidirectional  tensile  stress  (strain), 

'S  =  ultimate  unidirectional  compressive  stress  (strain), 
cr=  principal  or  local  stress  (strain). 

In  stress  (strain)  space,  a  general  maximum  stress  (strain)  failure  criteria  is  uncoupled, 
i.e.,  failure  is  unaffected  by  interaction  of  the  stresses;  see  Fig.  G.l. 
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Fig.  G.l  Stress  (strain)  space. 

In  actuality,  global  stresses  are  coupled  during  principal  stress  or  local  stress  calculations. 

V 

When  maximum  stresses  (strains)  are  examined  in  strain  (stress)  space,  the  strains 
(stresses)  are  coupled  by  Hook’s  law: 

sx  =-4(<x,  ~ua2  -uo-3) 

E 

s 2  =4(o-2-U<T1  -VC Tj) 

E 

*3  -^1  -^2) 

E 

G.3.  Ductile  Material  Failure 

Failure  of  ductile  materials  involves  predicting  transition  from  elastic  to  plastic  behavior 
and  the  concept  of  a  yield  surface.  Tresca  and  von  Mises  are  two  widely  used  ductile 

failure  criteria. 

G.3.1.  Maximum  Shear  Stress  Failure  Criterion 

Tresca,  or  maximum  shear  stress,  criterion  states  that  yielding  begins  in  a  material  point 
under  a  multi-axial  stress  state  when  maximum  shear  stress  equals  maximum  shear  stress 
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of  a  point  at  yield  in  a  uni-axial  tension  or  compression  test  of  the  same  material.  The 


material  fails  if: 


where 


>S 


s  > 


(G.7) 


and  w  is  the  largest  of 


^  _  ^uniaxial 

s  2 


(j .  =  principle  stress 


G.3.2.  Distortional  Energy  Density  (von  Mises) 

The  von  Mises  criterion,  also  known  as  distortional  energy  density  or  maximum 
octahedral  shear-stress  criterion,  uses  the  concept  that  strain  energy  at  a  point  can  be 
broken  into  two  components.  Part  of  the  energy  is  used  for  volumetric  change  and  the 
remainder  in  element  distortion.  The  von  Mises  criterion  states  that  yielding  begins  in  a 
material  point  under  a  multi-axial  stress  state  when  distortional  strain  energy  density 
equals  distortional  strain  energy  density  of  a  point  at  yield  in  a  uni-axial  tension  or 
compression  test  of  the  same  material.  The  von  Mises  criterion  can  be  stated  in  several 
forms.  Failure  occurs  when: 

Uo^U^  (G-8) 

where 

(<r, -a,f  +(<T,-<TJ)!+(<rs-<T,)’ 

U’~ - \2G 

Uoy  2G  * 

G  =  Shear  modulus 

also 
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1  3  2 

U  D  =  J  2  U  D  ~  T  octahedral 


with  the  octahedral  shear  stress  and  second  deviatoric  stress  invariant  defined  as 

2  _  (q-,  zEjf  tith  ~g3)2  +^3 

^ octahedral  C) 

(cr,  —  o';)2  +(g~2  —  <^3)  +(cr3~t7i) 

6 


^2  = 


Implied  in  both  von  Mises  and  Tresca  failure  criterion  is  that  hydrostatic  stress, 

=(°'i  +<t2  +o-3)/3’ 

has  no  effect  on  yielding. 

G.3.3.  General  Quadratic  Interaction  Criteria 

The  isotropic  von  Mises  failure  criterion  is  a  special  case  of  a  general  form  called 
Quadratic  Interaction  Criteria.  So  named  because  they  include  terms  to  account  for 
interaction  between  the  stress  components.  Hill  [5]  generalized  the  von  Mises  criterion  to 
include  the  effects  of  anisotropic  behavior.  The  Hill  criterion  is  given  by: 


where 


UD>  1 


UD  =  -cr2)2  +F2(o-2 -ct3)2  +  F3(ct3 -cr,)2 

+  ^23  (^23  +  ^32  )2  +Fn  (CT13  +  *ilY  +  FU  (ff12  +  ^21  Y 


(G.9) 


and 
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11  1  _  _J . J _ L_  9  F  =  —  +  — - — 

2  =  s222  «  sy  3  Si  s * 


2  ^=TT 

°i2 


2F13  =  -r- 

13  n2 
^13 


2^3=  TT 

*^23 


~  uniaxial  yield  strengths  in  the  1 ,2,3  directions, 
S  ~  Shear  yield  strengths  in  the  ij  direction. 


Stresses  are  referenced  to  the  principal  axes  of  anisotropy  (principal  material  directions), 
i.e.,  this  criterion  cannot  be  used  in  principal  (coordinate  invariant)  stress  space.  As 
presented,  Hill’s  criterion  assumes  no  difference  between  tensile  and  compressive  yield 
strengths. 


The  most  general  form  of  quadratic  interaction  criteria  [2]  as  suggested  by  Gol’denblat 
and  Kopnov  [18]  and  proposed  by  Tsai  and  Wu  [19]  is  of  the  form: 

.  Fn  +  FjrPjZ  1  <G-10) 

where  contracted  notation  is  used  (ij  =  1  to  6)  and  F,  and  Fy  are  experimentally 
determined  strength  tensors.  Hoffinan  [20]  has  suggested  keeping  the  linear  terms  in  the 
normal  stresses  as  a  way  to  account  for  differences  in  tensile  and  compressive  strengths. 

G.3 .4.  Strain  Energy  Density  Criterion 

In  elastic  non-linear  materials,  both  volumetric  and  distortional  changes  may  have  an 
effect  on  failure  hence  the  strain  energy  density  failure  criterion  may  be  appropriate.  The 
strain  energy  density  criterion  states  that  yielding  begins  in  a  material  point  under  a 
multi-axial  stress  state  when  total  strain  energy  density  equals  total  strain  energy  density 
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of  a  point  at  yield  in  a  uni-axial  tension  or  compression  test  of  the  same  material.  Failure 
under  the  strain  energy  density  criterion  occurs  when: 

i/.-u*  (GU) 

where 

U°  i8Jf  12G 

-_2 

rr  __  ^ uniaxial 

Uoy~  2E 
K  ~  Bulk  Modulus. 

Noting  that  the  first  term  in  the  expression  for  U0  is  strain  energy  density  associated  with 
pure  volume  change  and  the  second  term  is  distortional  strain  energy  density. 

G.4.  Discussion 

At  this  time  there  appears  to  be  little  to  no  advantage  in  using,  in  a  general  manner,  one 
failure  criteria  over  another.  The  accuracy  of  a  particular  criteria  is  likely  to  be  case 
dependent.  Arguments  for  such  a  statement  begin  with  an  observation  that  none  of  the 
available  failure  theories  have  been  shown  to  accurately  predict  failure  for  all  materials 
and  loading  conditions.  Furthermore,  there  is  no  consensus  within  the  engineering 
community  as  to  which  theory  is  best  [2]. 

Ultimately  the  goal  of  a  failure  prediction  methodology  is  to  extrapolate  results  of  an 
easy-to-conduct  laboratory  experiment,  simple  test  specimen  geometry  and  uniaxial 
loading,  to  real-world  applications  with  complex  structural  geometiy  and  multiaxial 
loading.  All  failure  criteria  mentioned  above  are  semi-empirical,  phenomenological 
attempts  to  develop  analytical  models  to  describe  experimental  observations. 
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Phenomenological,  because  only  the  occurrence  of  failure  is  predicted,  not  mode  of 
failure.  Failure  details  such  as  fiber  pullout,  fiber  breakage,  fiber  micro-buckling,  matrix 
cracking,  delamination,  etc.,  are  not  discerned.  Wu  [2,6]  states  that  a  large  experimental 
data  base  alone  could  form  the  basis  for  an  empirical  failure  criterion,  but  semi-empirical 
mathematical  models  are  preferable  because  they  can  reduce  the  number  of  required 
experiments  and  provide  a  more  systematic  approach  to  design  [2]. 
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